From 8de51122f13db7193a81318a53d1252d0a7e88a1 Mon Sep 17 00:00:00 2001 From: akif Date: Wed, 29 Jul 2026 18:23:53 +0200 Subject: [PATCH 01/20] feat(mip): bound presolve by rounds, badge size and work units Papilo and the probing cache were bounded only by a fraction of the wall clock (presolve_time_ratio * time_limit, and max_time_on_probing). On wide instances that fraction expires mid-round, so which reductions survive depends on how fast the machine is: bab2 at a 15s budget got 1.5s of Papilo and a reduced problem of 1.82M nonzeros, versus 1.17M when allowed to finish. Bound presolve by structure instead. Papilo takes a round cap and a ceiling on probing.minbadgesize; the probing cache takes a budget in work units, counted exactly from propagation iterations and probes and folded in at the step barrier in a fixed order so it does not depend on thread count. The timers now carry the time actually remaining in the solve, which keeps presolve from overrunning the whole budget without truncating it early. presolve_budget_policy.hpp maps structural features onto these budgets, with one policy per hypothesis so a sweep isolates the mapping itself; `legacy` reproduces the previous behaviour as a baseline. The `size` policy encodes what measurement showed: probing.minbadgesize, not the round count, drives Papilo's cost, and on wide problems a large badge buys almost nothing (capping it at 32 left square41, supportcase6, nw04 and rail507 bit-identical while cutting presolve 2-26x), whereas narrow problems reduce measurably worse when capped. Add PRESOLVE_BUDGET, PRESOLVE_PAPILO, PRESOLVE_PAPILO_REDUCED and PRESOLVE_PROBING log lines carrying the features that went in, the budgets that came out, and the realised cost, so budgets can be refit offline. A work budget is worth 0.3s on one instance and 16s on another, so the log records units_per_s to translate one back into the other. Signed-off-by: akif --- .../utils/presolve_budget_sweep.py | 310 ++++++++++++++++++ .../mathematical_optimization/constants.h | 11 + .../mip/heuristics_hyper_params.hpp | 23 +- cpp/src/math_optimization/solver_settings.cu | 7 + .../diversity/diversity_manager.cu | 41 ++- .../mip_heuristics/presolve/multi_probe.cu | 1 + .../mip_heuristics/presolve/multi_probe.cuh | 6 + .../presolve/presolve_budget_policy.hpp | 215 ++++++++++++ .../mip_heuristics/presolve/probing_cache.cu | 77 ++++- .../mip_heuristics/presolve/probing_cache.cuh | 33 +- .../presolve/third_party_presolve.cpp | 73 ++++- .../presolve/third_party_presolve.hpp | 12 +- cpp/src/mip_heuristics/solve.cu | 61 +++- 13 files changed, 824 insertions(+), 46 deletions(-) create mode 100644 benchmarks/linear_programming/utils/presolve_budget_sweep.py create mode 100644 cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp diff --git a/benchmarks/linear_programming/utils/presolve_budget_sweep.py b/benchmarks/linear_programming/utils/presolve_budget_sweep.py new file mode 100644 index 0000000000..0c5e82298a --- /dev/null +++ b/benchmarks/linear_programming/utils/presolve_budget_sweep.py @@ -0,0 +1,310 @@ +# SPDX-FileCopyrightText: Copyright (c) 2026, NVIDIA CORPORATION & AFFILIATES. All rights reserved. +# SPDX-License-Identifier: Apache-2.0 +r"""Sweep the MIP presolve budget policies over a set of instances. + +Runs cuopt_cli once per (instance, policy), scrapes the PRESOLVE_* log lines plus the final result +line, and writes one CSV row per run. The CSV carries both the structural features each budget was +derived from and what the budget actually spent, which is what a fit of +"budget <- problem dimensions and structure" needs. + +Sweeping whole policies:: + + python presolve_budget_sweep.py --cli ./cpp/build/cuopt_cli \ + --dataset-dir datasets/mip/miplib2017 --time-limit 300 --policies 0 2 \ + --out /tmp/presolve_sweep.csv --log-dir /tmp/presolve_sweep_logs + +Sweeping one knob instead gives the mapping curve from a wall-clock limit onto a +round / badge / work-unit budget:: + + python presolve_budget_sweep.py --grid-param mip_hyper_heuristic_presolve_max_rounds \ + --grid-values 1 2 3 5 10 20 30 50 --time-limit 2 --out /tmp/rounds_map.csv +""" + +import argparse +import csv +import os +import re +import subprocess +import sys +import time + +# Chosen to span the structural axes the policies key on: nnz, average row length, and binary +# fraction. bab2 / supportcase6 / square41 are the instances whose presolve ran away unbounded and +# motivated the budgets in the first place, so they are the regression cases. +DEFAULT_INSTANCES = [ + "bab2", + "supportcase6", + "square41", + "air05", + "30n20b8", + "gen-ip054", + "nw04", + "rail507", + "seymour", + "mzzv11", + "roll3000", + "ns1208400", + "glass4", + "timtab1", + "enlight_hard", + "sp97ar", +] + +POLICY_NAMES = { + 0: "legacy", + 1: "fixed", + 2: "size", + 3: "density", + 4: "binary", + 5: "combined", + 6: "manual", +} + +KV_RE = re.compile(r"(\w+)=(-?[\w.+-]+)") +EXPLORED_RE = re.compile( + r"Explored (\d+) nodes \((\d+) simplex iterations\) in ([\d.]+)s" +) +OBJ_RE = re.compile( + r"Best objective ([-\d.eE+]+), best bound ([-\d.eE+]+), gap ([-\d.eE+]+|inf)%" +) + +STATUS_MARKERS = [ + ("Optimal solution found", "Optimal"), + ("Time limit reached", "TimeLimit"), + ("Work limit reached", "WorkLimit"), + ("Problem is infeasible", "Infeasible"), + ("Problem is unbounded", "Unbounded"), + ("No solution found", "NoSolution"), + # A claim of integer infeasibility on an instance known to be feasible means a reduction was + # unsound, so it must be distinguishable from simply not finding a solution in time. + ("Problem has no integer feasible solution", "NoIntegerFeasible"), +] + + +def parse_kv(line, prefix): + """Pull every key=value token that follows `prefix` on `line`.""" + idx = line.find(prefix) + if idx < 0: + return {} + return dict(KV_RE.findall(line[idx + len(prefix) :])) + + +def parse_log(text): + row = {} + for line in text.splitlines(): + if "PRESOLVE_BUDGET stage=PAPILO" in line: + for k, v in parse_kv(line, "PRESOLVE_BUDGET").items(): + row[f"papilo_{k}"] = v + elif "PRESOLVE_BUDGET stage=PROBING" in line: + for k, v in parse_kv(line, "PRESOLVE_BUDGET").items(): + row[f"probing_{k}"] = v + elif "PRESOLVE_PAPILO_REDUCED" in line: + for k, v in parse_kv(line, "PRESOLVE_PAPILO_REDUCED").items(): + row[f"reduced_{k}"] = v + elif "PRESOLVE_PAPILO wall=" in line: + for k, v in parse_kv(line, "PRESOLVE_PAPILO").items(): + row["papilo_wall" if k == "wall" else f"papilo_{k}"] = v + elif "PRESOLVE_PROBING_WALL" in line: + row["probing_wall"] = parse_kv(line, "PRESOLVE_PROBING_WALL").get( + "wall" + ) + elif "PRESOLVE_PROBING probes=" in line: + for k, v in parse_kv(line, "PRESOLVE_PROBING").items(): + row[f"spent_{k}"] = v + elif "Probing-cache step disabled" in line: + row["probing_disabled"] = "1" + + m = EXPLORED_RE.search(line) + if m: + row["nodes"], row["simplex_iters"], row["solve_wall"] = m.groups() + m = OBJ_RE.search(line) + if m: + row["objective"], row["bound"], row["gap_pct"] = m.groups() + for marker, status in STATUS_MARKERS: + if marker in line: + row["status"] = status + return row + + +def as_text(stream): + """Decode a stream, since subprocess hands back bytes on the timeout path even with text=True.""" + if stream is None: + return "" + if isinstance(stream, bytes): + return stream.decode("utf-8", "replace") + return stream + + +def run_one(args, instance, policy, config_path, grid_value=None): + with open(config_path, "w") as fh: + fh.write(f"mip_hyper_heuristic_presolve_budget_policy = {policy}\n") + if grid_value is not None: + fh.write(f"{args.grid_param} = {grid_value}\n") + for extra in args.param: + fh.write(extra.replace(":", " = ", 1) + "\n") + + cmd = [ + args.cli, + os.path.join(args.dataset_dir, instance + ".mps"), + "--time-limit", + str(args.time_limit), + "--params-file", + config_path, + ] + if args.determinism: + cmd += ["--mip-determinism-mode", "1"] + + t0 = time.time() + timed_out = False + try: + proc = subprocess.run( + cmd, + capture_output=True, + text=True, + timeout=args.timeout, + check=False, + ) + output = proc.stdout + proc.stderr + returncode = proc.returncode + except subprocess.TimeoutExpired as exc: + output = as_text(exc.stdout) + as_text(exc.stderr) + returncode = -1 + timed_out = True + wall = time.time() - t0 + + row = parse_log(output) + row.update( + instance=instance, + policy=policy, + policy_name=POLICY_NAMES.get(policy, str(policy)), + harness_wall=f"{wall:.2f}", + harness_timeout=int(timed_out), + returncode=returncode, + ) + if grid_value is not None: + row["grid_param"] = args.grid_param + row["grid_value"] = grid_value + if timed_out: + row.setdefault("status", "HarnessTimeout") + + if args.log_dir: + os.makedirs(args.log_dir, exist_ok=True) + with open( + os.path.join(args.log_dir, f"{instance}.p{policy}.log"), "w" + ) as fh: + fh.write(output) + return row + + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--cli", default="./cpp/build/cuopt_cli") + ap.add_argument("--dataset-dir", default="datasets/mip/miplib2017") + ap.add_argument("--instances", nargs="*", default=DEFAULT_INSTANCES) + ap.add_argument( + "--policies", nargs="*", type=int, default=[0, 1, 2, 3, 4, 5] + ) + ap.add_argument("--time-limit", type=float, default=30.0) + ap.add_argument( + "--timeout", + type=float, + default=300.0, + help="hard wall cap per run; the legacy policy leaves presolve unbounded", + ) + ap.add_argument("--determinism", action="store_true") + ap.add_argument( + "--param", + action="append", + default=[], + help="extra config entry as key:value, repeatable", + ) + ap.add_argument( + "--grid-param", + default="", + help="sweep this single hyper-parameter under the manual policy instead of sweeping " + "policies; this is what yields the wall-limit -> rounds / work-unit mapping curve", + ) + ap.add_argument("--grid-values", nargs="*", default=[]) + ap.add_argument("--out", default="presolve_sweep.csv") + ap.add_argument("--log-dir", default="") + args = ap.parse_args() + + config_path = "/tmp/presolve_budget_sweep.config" + rows = [] + # A grid sweeps one knob under the manual policy; otherwise the variant axis is the policy. + if args.grid_param: + variants = [(6, v) for v in args.grid_values] + else: + variants = [(p, None) for p in args.policies] + + total = len(args.instances) * len(variants) + done = 0 + for instance in args.instances: + path = os.path.join(args.dataset_dir, instance + ".mps") + if not os.path.exists(path): + print(f"SKIP missing {path}", flush=True) + continue + for policy, grid_value in variants: + done += 1 + # A sweep is long enough that losing all of it to one bad run is the worst outcome. + try: + row = run_one(args, instance, policy, config_path, grid_value) + except Exception as exc: # noqa: BLE001 + row = { + "instance": instance, + "policy": policy, + "policy_name": POLICY_NAMES.get(policy, str(policy)), + "grid_value": grid_value, + "status": "HarnessError", + "harness_wall": "0", + "harness_error": repr(exc), + } + rows.append(row) + variant = ( + f"{args.grid_param.split('_')[-1]}={grid_value}" + if grid_value is not None + else f"p{policy} {row.get('policy_name', '')}" + ) + print( + f"[{done}/{total}] {instance:16s} {variant:14s}" + f" status={row.get('status', '?'):18s}" + f" papilo_wall={row.get('papilo_wall', '-'):>8s}" + f" red_vars={row.get('reduced_nvars', '-'):>8s}" + f" red_nnz={row.get('reduced_nnz', '-'):>9s}" + f" probing_wall={row.get('probing_wall', '-'):>8s}" + f" probes={row.get('spent_probes', '-'):>7s}" + f" work={row.get('spent_work', '-'):>9s}", + flush=True, + ) + # Written incrementally so a long sweep is inspectable while it runs. + write_csv(args.out, rows) + + write_csv(args.out, rows) + print(f"\nwrote {len(rows)} rows to {args.out}") + + +def write_csv(path, rows): + if not rows: + return + fields = [] + for row in rows: + for key in row: + if key not in fields: + fields.append(key) + lead = [ + "instance", + "policy", + "policy_name", + "grid_param", + "grid_value", + "status", + ] + fields = lead + [f for f in fields if f not in lead] + with open(path, "w", newline="") as fh: + writer = csv.DictWriter(fh, fieldnames=fields, extrasaction="ignore") + writer.writeheader() + writer.writerows(rows) + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/cpp/include/cuopt/mathematical_optimization/constants.h b/cpp/include/cuopt/mathematical_optimization/constants.h index 9592389dea..da92ac8733 100644 --- a/cpp/include/cuopt/mathematical_optimization/constants.h +++ b/cpp/include/cuopt/mathematical_optimization/constants.h @@ -99,6 +99,17 @@ #define CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS "mip_hyper_heuristic_num_cpufj_threads" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO "mip_hyper_heuristic_presolve_time_ratio" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME "mip_hyper_heuristic_presolve_max_time" +#define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY \ + "mip_hyper_heuristic_presolve_budget_policy" +#define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS "mip_hyper_heuristic_presolve_max_rounds" +#define CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE \ + "mip_hyper_heuristic_papilo_probing_max_badgesize" +#define CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT \ + "mip_hyper_heuristic_cuopt_presolve_work_limit" +#define CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE "mip_hyper_heuristic_probing_step_size" +#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK \ + "mip_hyper_heuristic_probe_host_overhead_work" +#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK "mip_hyper_heuristic_probe_iter_work" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO "mip_hyper_heuristic_root_lp_time_ratio" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME "mip_hyper_heuristic_root_lp_max_time" #define CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT "mip_hyper_heuristic_rins_time_limit" diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index e3fd891e83..e2ed78e2e4 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -19,10 +19,25 @@ namespace cuopt::mathematical_optimization { */ template struct mip_heuristics_hyper_params_t { - i_t population_size = 32; // max solutions in pool - i_t num_cpufj_threads = 8; // parallel CPU FJ climbers - f_t presolve_time_ratio = 0.1; // fraction of total time for presolve - f_t presolve_max_time = 60.0; // hard cap on presolve seconds + i_t population_size = 32; // max solutions in pool + i_t num_cpufj_threads = 8; // parallel CPU FJ climbers + f_t presolve_time_ratio = 0.1; // fraction of total time for presolve + f_t presolve_max_time = 60.0; // hard cap on presolve seconds + + // Presolve budgeting. presolve_budget_policy selects how the four knobs below are derived from + // the problem's dimensions and structure (see presolve_budget_policy.hpp); the values here are + // the defaults the policy starts from and the literal values used by the `manual` policy. + i_t presolve_budget_policy = 1; // presolve_budget_policy_t + i_t presolve_max_rounds = 30; // Papilo presolve rounds cap (<=0 = Papilo default) + i_t papilo_probing_max_badgesize = 1024; // ceiling on Papilo's probing.minbadgesize + f_t cuopt_presolve_work_limit = 30.0; // probing-cache budget, work units + i_t probing_step_size = 512; // probed vars between work-budget checks + // Weights of the probing-cache work model. Work units measure probing effort reproducibly; they + // are not an estimate of elapsed time, so the effort-per-second they correspond to legitimately + // differs between instances. + f_t probe_host_overhead_work = 0.02; // charged per probed variable + f_t probe_iter_work = 0.01; // charged per multi-probe propagation iteration + f_t root_lp_time_ratio = 0.1; // fraction of total time for root LP f_t root_lp_max_time = 15.0; // hard cap on root LP seconds f_t rins_time_limit = 3.0; // per-call RINS sub-MIP time diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index 666a12bb91..bedc60a053 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -106,6 +106,9 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO, &mip_settings.heuristic_params.presolve_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for presolve"}, {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME, &mip_settings.heuristic_params.presolve_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(60.0), "hard cap on presolve seconds"}, + {CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT, &mip_settings.heuristic_params.cuopt_presolve_work_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(30.0), "probing-cache budget in work units (manual policy, or ceiling for derived policies)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK, &mip_settings.heuristic_params.probe_host_overhead_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.02), "work units charged per probed variable (host overhead)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK, &mip_settings.heuristic_params.probe_iter_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.01), "work units charged per multi-probe propagation iteration"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, @@ -167,6 +170,10 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") {CUOPT_MIP_HYPER_HEURISTIC_POPULATION_SIZE, &mip_settings.heuristic_params.population_size, 1, std::numeric_limits::max(), 32, "max solutions in pool"}, {CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS, &mip_settings.heuristic_params.num_cpufj_threads, 0, std::numeric_limits::max(), 8, "parallel CPU FJ climbers"}, + {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY, &mip_settings.heuristic_params.presolve_budget_policy, 0, 6, 1, "how presolve budgets are derived: 0=legacy 1=fixed 2=size 3=density 4=binary 5=combined 6=manual"}, + {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS, &mip_settings.heuristic_params.presolve_max_rounds, -1, std::numeric_limits::max(), 30, "Papilo presolve rounds cap (<=0 keeps Papilo default)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE, &mip_settings.heuristic_params.papilo_probing_max_badgesize, -1, std::numeric_limits::max(), 1024, "ceiling on Papilo probing.minbadgesize (<=0 leaves it uncapped)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE, &mip_settings.heuristic_params.probing_step_size, 1, std::numeric_limits::max(), 512, "probed variables between probing-cache work-budget checks"}, {CUOPT_MIP_HYPER_HEURISTIC_STAGNATION_TRIGGER, &mip_settings.heuristic_params.stagnation_trigger, 1, std::numeric_limits::max(), 3, "FP loops w/o improvement before recombination"}, {CUOPT_MIP_HYPER_HEURISTIC_MAX_ITERS_WITHOUT_IMPROVEMENT, &mip_settings.heuristic_params.max_iterations_without_improvement, 1, std::numeric_limits::max(), 8, "diversity step depth after stagnation"}, {CUOPT_MIP_HYPER_HEURISTIC_N_OF_MINIMUMS_FOR_EXIT, &mip_settings.heuristic_params.n_of_minimums_for_exit, 1, std::numeric_limits::max(), 7000, "FJ baseline local-minima exit threshold"}, diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 61c90944f4..02b3fc4c62 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -20,6 +20,9 @@ #include +#include +#include +#include #include #include @@ -294,10 +297,19 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ if (termination_criterion_t::NO_UPDATE != term_crit) { ls.constraint_prop.bounds_update.set_updated_bounds(*problem_ptr); } - bool run_probing_cache = !fj_only_run; - // Don't run probing cache in deterministic mode yet as neither B&B nor CPUFJ need it - // and it doesn't make use of work units yet - if (context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC) { run_probing_cache = false; } + const auto& hp = context.settings.heuristic_params; + const bool deterministic = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC; + const auto probing_features = probing_presolve_features(*problem_ptr); + const auto probing_budget = evaluate_presolve_budget(hp, probing_features); + bool run_probing_cache = !fj_only_run; + // Under the legacy policy the probing cache carries no work budget, so in deterministic mode + // nothing would bound it: the wall clock is infinite there. Keep it off, which is also what makes + // the legacy policy an unchanged baseline. + if (deterministic && + probing_budget.probing_work_limit == std::numeric_limits::infinity()) { + CUOPT_LOG_INFO("Probing-cache step disabled: deterministic mode with no work budget"); + run_probing_cache = false; + } // Allow the user to disable the probing-cache step of cuOpt's internal presolve // independently of the higher-level presolver setting. if (!context.settings.probing) { @@ -305,13 +317,22 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ run_probing_cache = false; } if (run_probing_cache) { - // Run probing cache before trivial presolve to discover variable implications - const f_t max_time_on_probing = diversity_config.max_time_on_probing; - f_t time_for_probing_cache = std::min(max_time_on_probing, time_limit); - timer_t probing_timer{time_for_probing_cache}; + log_presolve_budget("PROBING", hp.presolve_budget_policy, probing_features, probing_budget); + // Run probing cache before trivial presolve to discover variable implications. The work budget + // decides how much probing happens; the timer only carries the time actually left in the solve, + // so it cannot overrun the whole budget but also is not what shapes the amount of probing. + timer_t probing_timer{global_timer.remaining_time()}; + const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem - bool problem_is_infeasible = - compute_probing_cache(ls.constraint_prop.bounds_update, *problem_ptr, probing_timer); + bool problem_is_infeasible = compute_probing_cache(ls.constraint_prop.bounds_update, + *problem_ptr, + probing_timer, + probing_budget.probing_work_limit, + (size_t)probing_budget.probing_step_size); + problem_ptr->handle_ptr->sync_stream(); + CUOPT_LOG_INFO( + "PRESOLVE_PROBING_WALL wall=%.3f", + std::chrono::duration(std::chrono::steady_clock::now() - probing_t0).count()); if (problem_is_infeasible) { return false; } } const bool remap_cache_ids = true; diff --git a/cpp/src/mip_heuristics/presolve/multi_probe.cu b/cpp/src/mip_heuristics/presolve/multi_probe.cu index f1adf28650..394d89f580 100644 --- a/cpp/src/mip_heuristics/presolve/multi_probe.cu +++ b/cpp/src/mip_heuristics/presolve/multi_probe.cu @@ -306,6 +306,7 @@ termination_criterion_t multi_probe_t::bound_update_loop(problem_tsync_stream(); if (compute_stats) { upd_0.init_changed_constraints(handle_ptr); diff --git a/cpp/src/mip_heuristics/presolve/multi_probe.cuh b/cpp/src/mip_heuristics/presolve/multi_probe.cuh index b4227281b6..2f2f2f038d 100644 --- a/cpp/src/mip_heuristics/presolve/multi_probe.cuh +++ b/cpp/src/mip_heuristics/presolve/multi_probe.cuh @@ -75,6 +75,12 @@ class multi_probe_t { bool skip_0; bool skip_1; settings_t settings; + // When set, the number of propagation iterations run is accumulated here rather than on any + // shared counter. The probing cache runs one multi_probe_t per OMP task, so a shared counter + // would be an unsynchronized read-modify-write and would make the resulting budget + // nondeterministic. The owner folds the per-task counts in at a barrier, in a fixed order, and + // applies its own cost model. + double* local_iter_accumulator = nullptr; bool compute_stats = true; bool init_changed_constraints = true; i_t infeas_constraints_count_0 = 0; diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp new file mode 100644 index 0000000000..891b598cb2 --- /dev/null +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -0,0 +1,215 @@ +/* clang-format off */ +/* + * SPDX-FileCopyrightText: Copyright (c) 2026, NVIDIA CORPORATION & AFFILIATES. All rights reserved. + * SPDX-License-Identifier: Apache-2.0 + */ +/* clang-format on */ + +#pragma once + +#include + +#include + +#include +#include +#include + +namespace cuopt::mathematical_optimization::mip { + +// Competing hypotheses for how much presolve effort a problem deserves. Each policy maps the +// structural features below onto the same four budgets, so a sweep over policies on a fixed +// instance set isolates the effect of the mapping itself. `legacy` reproduces pre-budget behaviour +// and is the baseline to measure against; `manual` reads the hyper-parameters verbatim so a +// specific point can be pinned from the command line. +enum class presolve_budget_policy_t : int { + legacy = 0, + fixed = 1, + size = 2, + density = 3, + binary = 4, + combined = 5, + manual = 6, +}; + +inline const char* presolve_budget_policy_name(int policy) +{ + switch (static_cast(policy)) { + case presolve_budget_policy_t::legacy: return "legacy"; + case presolve_budget_policy_t::fixed: return "fixed"; + case presolve_budget_policy_t::size: return "size"; + case presolve_budget_policy_t::density: return "density"; + case presolve_budget_policy_t::binary: return "binary"; + case presolve_budget_policy_t::combined: return "combined"; + case presolve_budget_policy_t::manual: return "manual"; + default: return "unknown"; + } +} + +// Dimensions plus cheap structural ratios. Both presolve stages populate this from whatever problem +// representation they hold: Papilo from the original problem before any reduction, the cuOpt +// probing cache from the Papilo-reduced problem. The two therefore see different feature values for +// the same instance, which is intended -- each budget should follow the problem it actually +// operates on. +struct presolve_features_t { + double n_vars{0}; + double n_cons{0}; + double nnz{0}; + double n_int{0}; + double n_bin{0}; + double max_row_len{0}; + + double avg_row_len() const { return n_cons > 0 ? nnz / n_cons : 0.0; } + double avg_col_len() const { return n_vars > 0 ? nnz / n_vars : 0.0; } + double density() const { return (n_vars > 0 && n_cons > 0) ? nnz / (n_vars * n_cons) : 0.0; } + double int_frac() const { return n_vars > 0 ? n_int / n_vars : 0.0; } + double bin_frac() const { return n_vars > 0 ? n_bin / n_vars : 0.0; } +}; + +struct presolve_budget_t { + // <=0 leaves Papilo's own default (unlimited rounds). + int papilo_max_rounds{-1}; + // <=0 leaves probing.minbadgesize uncapped at max(ncols/2, 32). + int papilo_max_badgesize{-1}; + // Probing-cache budget in work units: a reproducible count of probing effort, not a time + // estimate. + double probing_work_limit{std::numeric_limits::infinity()}; + // Probed variables per step, i.e. the granularity at which the budget can be enforced. + int probing_step_size{2048}; +}; + +namespace detail { + +inline double clamp_d(double v, double lo, double hi) { return std::min(std::max(v, lo), hi); } + +inline int clamp_i(double v, int lo, int hi) +{ + return static_cast(clamp_d(std::round(v), lo, hi)); +} + +} // namespace detail + +template +presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t& hp, + const presolve_features_t& feat) +{ + using detail::clamp_d; + using detail::clamp_i; + + presolve_budget_t b{}; + + const double nnz = std::max(feat.nnz, 1.0); + const double arl = std::max(feat.avg_row_len(), 1.0); + const double n_int = std::max(feat.n_int, 1.0); + const double n_bin = std::max(feat.n_bin, 1.0); + const double bf = feat.bin_frac(); + + // Lets a single knob scale every derived policy without recompiling: 1.0 at the default of 30. + const double work_scale = static_cast(hp.cuopt_presolve_work_limit) / 30.0; + double raw_work = 30.0; + + switch (static_cast(hp.presolve_budget_policy)) { + case presolve_budget_policy_t::legacy: + b.papilo_max_rounds = -1; + b.papilo_max_badgesize = -1; + b.probing_work_limit = std::numeric_limits::infinity(); + b.probing_step_size = 2048; + return b; + + case presolve_budget_policy_t::fixed: + b.papilo_max_rounds = 30; + b.papilo_max_badgesize = 1024; + b.probing_work_limit = 30.0; + b.probing_step_size = 512; + return b; + + case presolve_budget_policy_t::manual: + b.papilo_max_rounds = hp.presolve_max_rounds; + b.papilo_max_badgesize = hp.papilo_probing_max_badgesize; + b.probing_work_limit = hp.cuopt_presolve_work_limit; + b.probing_step_size = hp.probing_step_size; + return b; + + // Measured rule. probing.minbadgesize, not the round count, is what drives Papilo's cost, and + // on wide problems a large badge buys almost nothing: capping it at 32 left the reduced problem + // bit-identical on square41/supportcase6/nw04/rail507 while cutting presolve 2-26x. On narrow + // problems the opposite holds (mzzv11, 30n20b8, air05 all reduce measurably worse at 32), so + // the cap is applied by width only. Rounds are kept non-binding on wide problems -- everything + // measured saturates well before 50 -- and capped on narrow ones, where mzzv11 keeps growing. + case presolve_budget_policy_t::size: { + const bool wide = feat.n_vars > 2.0e4; + b.papilo_max_rounds = wide ? 50 : 20; + b.papilo_max_badgesize = wide ? 32 : -1; + raw_work = 120.0; + b.probing_step_size = 512; + break; + } + + // A single propagation sweep costs roughly one pass over each touched row, so long rows make + // every probe more expensive: buy fewer of them, and size the badge so one badge's working + // limit (~2*nnz in Papilo) stays roughly constant. + case presolve_budget_policy_t::density: + b.papilo_max_rounds = arl <= 10 ? 40 : arl <= 50 ? 20 : 8; + b.papilo_max_badgesize = clamp_i(2.0e6 / arl, 32, 4096); + raw_work = 30.0 * (10.0 / arl); + b.probing_step_size = arl <= 50 ? 1024 : 256; + break; + + // Probing and clique merging only pay off on binaries, so scale with how many there are and how + // much of the problem they make up. + case presolve_budget_policy_t::binary: + b.papilo_max_rounds = bf >= 0.9 ? 50 : bf >= 0.5 ? 30 : 15; + b.papilo_max_badgesize = clamp_i(std::max(n_bin / 2.0, 32.0), 32, 1024); + raw_work = (2.0 / 3.0) * std::sqrt(n_bin); + b.probing_step_size = 512; + break; + + // Multiplicative over the three effects above, so no single feature can dominate the budget. + case presolve_budget_policy_t::combined: { + const double size_f = clamp_d(1.0e5 / nnz, 0.25, 2.0); + const double density_f = clamp_d(10.0 / arl, 0.25, 2.0); + const double binary_f = clamp_d(0.5 + bf, 0.5, 1.5); + const double factor = size_f * density_f * binary_f; + b.papilo_max_rounds = clamp_i(30.0 * factor, 5, 60); + b.papilo_max_badgesize = clamp_i(2.0e6 / arl, 32, 2048); + raw_work = 30.0 * factor; + b.probing_step_size = clamp_i(512.0 * density_f, 128, 2048); + break; + } + } + + b.probing_work_limit = work_scale * clamp_d(raw_work, 5.0, 120.0); + return b; +} + +// One line per presolve stage carrying the features that went in and the budgets that came out, so +// a sweep can be regressed offline without re-deriving anything from the solver. +inline void log_presolve_budget(const char* stage, + int policy, + const presolve_features_t& f, + const presolve_budget_t& b) +{ + CUOPT_LOG_INFO( + "PRESOLVE_BUDGET stage=%s policy=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f nbin=%.0f " + "arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " + "rounds=%d badge=%d work=%.3f step=%d", + stage, + presolve_budget_policy_name(policy), + f.n_vars, + f.n_cons, + f.nnz, + f.n_int, + f.n_bin, + f.avg_row_len(), + f.avg_col_len(), + f.max_row_len, + f.density(), + f.int_frac(), + f.bin_frac(), + b.papilo_max_rounds, + b.papilo_max_badgesize, + b.probing_work_limit, + b.probing_step_size); +} + +} // namespace cuopt::mathematical_optimization::mip diff --git a/cpp/src/mip_heuristics/presolve/probing_cache.cu b/cpp/src/mip_heuristics/presolve/probing_cache.cu index fd4790479b..d3b7052600 100644 --- a/cpp/src/mip_heuristics/presolve/probing_cache.cu +++ b/cpp/src/mip_heuristics/presolve/probing_cache.cu @@ -21,6 +21,7 @@ #include #include +#include #include #include @@ -847,7 +848,9 @@ std::vector compute_priority_indices_by_implied_integers(problem_t bool compute_probing_cache(bound_presolve_t& bound_presolve, problem_t& problem, - timer_t timer) + timer_t timer, + double work_limit, + size_t step_size_hint) { raft::common::nvtx::range fun_scope("compute_probing_cache"); // we dont want to compute the probing cache for all variables for time and computation resources @@ -869,11 +872,16 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, std::vector>> modification_vector_pool(num_tasks); std::vector>> substitution_vector_pool(num_tasks); + // Each task counts its own propagation iterations; the counts are folded in at the step barrier + // below, where only one thread is running. + std::vector iter_accum_pool(num_tasks, 0.0); + // Initialize multi_probe_presolve_pool for (size_t i = 0; i < num_tasks; i++) { multi_probe_presolve_pool.emplace_back(bound_presolve.context); multi_probe_presolve_pool[i].resize(problem); - multi_probe_presolve_pool[i].compute_stats = true; + multi_probe_presolve_pool[i].compute_stats = true; + multi_probe_presolve_pool[i].local_iter_accumulator = &iter_accum_pool[i]; } // Atomic variables for tracking progress @@ -882,18 +890,39 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, std::atomic problem_is_infeasible(false); size_t last_it_implied_singletons = 0; bool early_exit = false; - const size_t step_size = min((size_t)2048, priority_indices.size()); + // Two additive terms with independently tunable weights, so a sweep can shift the budget between + // counting propagation iterations and counting probes. Both are exact counts, which is what keeps + // the budget reproducible. + const double iter_cost = (double)bound_presolve.context.settings.heuristic_params.probe_iter_work; + const double probe_cost = + (double)bound_presolve.context.settings.heuristic_params.probe_host_overhead_work; + // Only for the diagnostic below: a work budget buys wildly different amounts of time per + // instance, so the realised rate has to be recorded to translate a budget back into seconds + // afterwards. + const auto probing_t0 = std::chrono::steady_clock::now(); + double iters_done = 0.0; + size_t probes_done = 0; + double work_used = 0.0; + // Work is only folded in at the step barrier, so the step size is also the granularity at which + // the budget can be enforced: too large and a single step runs effectively unbudgeted. + const size_t step_size = min(step_size_hint, priority_indices.size()); // The pool buffers above were allocated on the main stream. // Each OMP thread below uses its own stream, so we must ensure all allocations // are visible before any per-thread kernel can reference that memory. problem.handle_ptr->sync_stream(); - CUOPT_LOG_INFO("Running probing cache with %zu tasks", num_tasks); + CUOPT_LOG_INFO( + "Running probing cache with %zu tasks (%zu candidate vars, work limit %.3f, step %zu)", + num_tasks, + priority_indices.size(), + work_limit, + step_size); // Main parallel loop for (size_t step_start = 0; step_start < priority_indices.size(); step_start += step_size) { if (timer.check_time_limit() || early_exit || problem_is_infeasible.load()) { break; } + if (work_used >= work_limit) { break; } size_t step_end = std::min(step_start + step_size, priority_indices.size()); #pragma omp taskloop num_tasks(num_tasks) default(shared) priority(CUOPT_DEFAULT_TASK_PRIORITY) @@ -926,8 +955,15 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, } } // implicit barrier that waits for all iterations to finish before proceeding - // TODO when we have determinism, check current threads work/time counter and filter queue - // items that are smaller or equal to that + // Single-threaded from here to the end of the step, so folding the per-task counts in a fixed + // order gives the same work_used for any thread count. + for (size_t t = 0; t < num_tasks; ++t) { + iters_done += iter_accum_pool[t]; + iter_accum_pool[t] = 0.0; + } + probes_done += step_end - step_start; + work_used = iters_done * iter_cost + (double)probes_done * probe_cost; + apply_modification_queue_to_problem(modification_vector_pool, problem); // copy host bounds again, because we changed some problem bounds raft::copy(h_var_bounds.data(), @@ -943,9 +979,28 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, } // end of step apply_substitution_queue_to_problem(substitution_vector_pool, problem); - CUOPT_LOG_DEBUG("Total number of cached probings %lu number of implied singletons %lu", - n_of_cached_probings.load(), - n_of_implied_singletons.load()); + const double probing_wall = + std::chrono::duration(std::chrono::steady_clock::now() - probing_t0).count(); + CUOPT_LOG_INFO( + "PRESOLVE_PROBING probes=%zu candidates=%zu iters=%.0f work=%.3f work_limit=%.3f step=%zu " + "iter_cost=%.5f probe_cost=%.5f wall=%.3f wall_limit=%.3f units_per_s=%.1f " + "budget_exhausted=%d early_exit=%d timed_out=%d cached=%lu implied_singletons=%lu", + probes_done, + priority_indices.size(), + iters_done, + work_used, + work_limit, + step_size, + iter_cost, + probe_cost, + probing_wall, + timer.get_time_limit(), + probing_wall > 0.0 ? work_used / probing_wall : 0.0, + (int)(work_used >= work_limit), + (int)early_exit, + (int)timer.check_time_limit(), + n_of_cached_probings.load(), + n_of_implied_singletons.load()); // restore the settings bound_presolve.settings = {}; return problem_is_infeasible.load(); @@ -954,7 +1009,9 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, #define INSTANTIATE(F_TYPE) \ template bool compute_probing_cache(bound_presolve_t & bound_presolve, \ problem_t & problem, \ - timer_t timer); \ + timer_t timer, \ + double work_limit, \ + size_t step_size_hint); \ template class probing_cache_t; #if MIP_INSTANTIATE_FLOAT diff --git a/cpp/src/mip_heuristics/presolve/probing_cache.cuh b/cpp/src/mip_heuristics/presolve/probing_cache.cuh index ec532febb9..74526d1ba2 100644 --- a/cpp/src/mip_heuristics/presolve/probing_cache.cuh +++ b/cpp/src/mip_heuristics/presolve/probing_cache.cuh @@ -11,8 +11,15 @@ #include +#include + +#include #include +#include +#include +#include + namespace cuopt::mathematical_optimization::mip { template @@ -116,9 +123,33 @@ class lb_probing_cache_t { std::unordered_map, 2>> probing_cache; }; +// Features of the Papilo-reduced problem, which is what the probing cache actually runs on. These +// differ from the features Papilo's own budget was derived from, sometimes by a lot. +template +presolve_features_t probing_presolve_features(problem_t const& problem) +{ + presolve_features_t f{}; + f.n_vars = problem.n_variables; + f.n_cons = problem.n_constraints; + f.nnz = problem.nnz; + f.n_int = problem.n_integer_vars; + f.n_bin = problem.n_binary_vars; + + auto h_offsets = cuopt::host_copy(problem.offsets, problem.handle_ptr->get_stream()); + for (size_t i = 0; i + 1 < h_offsets.size(); ++i) { + f.max_row_len = std::max(f.max_row_len, h_offsets[i + 1] - h_offsets[i]); + } + return f; +} + +// `work_limit` bounds probing in work units, checked at every step barrier and therefore +// independent of thread count and wall clock. `step_size_hint` is the number of variables probed +// per step, i.e. the granularity at which the budget can be enforced. template bool compute_probing_cache(bound_presolve_t& bound_presolve, problem_t& problem, - timer_t timer); + timer_t timer, + double work_limit = std::numeric_limits::infinity(), + size_t step_size_hint = 2048); } // namespace cuopt::mathematical_optimization::mip diff --git a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp index dc8410e7a1..614a23eb75 100644 --- a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp +++ b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp @@ -49,6 +49,7 @@ #include #include +#include #include #include #include @@ -709,11 +710,15 @@ void set_presolve_options(papilo::Presolve& presolver, f_t relative_tolerance, f_t time_limit, bool dual_postsolve, - i_t num_cpu_threads) + i_t num_cpu_threads, + i_t max_rounds) { presolver.getPresolveOptions().tlim = time_limit; presolver.getPresolveOptions().threads = num_cpu_threads; // user setting or 0 (automatic) presolver.getPresolveOptions().feastol = 1e-5; + // A round cap bounds presolve independently of the clock, which is the only thing that bounds it + // in deterministic mode where tlim is infinite. <=0 keeps Papilo's default of unlimited rounds. + if (max_rounds > 0) { presolver.getPresolveOptions().maxrounds = max_rounds; } if (dual_postsolve) { presolver.getPresolveOptions().componentsmaxint = -1; presolver.getPresolveOptions().detectlindep = 0; @@ -724,7 +729,8 @@ template void set_presolve_parameters(papilo::Presolve& presolver, problem_category_t category, int nrows, - int ncols) + int ncols, + int max_badgesize) { // It looks like a copy. But this copy has the pointers to relevant variables in papilo auto params = presolver.getParameters(); @@ -732,8 +738,13 @@ void set_presolve_parameters(papilo::Presolve& presolver, // Papilo has work unit measurements for probing. Because of this when the first batch fails to // produce any reductions, the algorithm stops. To avoid stopping the algorithm, we set a // minimum badge size to a huge value. The time limit makes sure that we exit if it takes too - // long + // long. + // An uncapped ncols/2 forces one probing pass to span the whole problem, so probing never + // reaches its work-based stop and runs unbounded on large MIPs whenever the clock is infinite. + // Capping the badge keeps it large enough to still find reductions while Papilo's per-badge + // working limit (~2*nnz) bounds a single pass. <=0 restores the uncapped behaviour. int min_badgesize = std::max(ncols / 2, 32); + if (max_badgesize > 0) { min_badgesize = std::min(min_badgesize, max_badgesize); } params.setParameter("probing.minbadgesize", min_badgesize); params.setParameter("cliquemerging.enabled", true); params.setParameter("cliquemerging.maxcalls", 50); @@ -818,7 +829,9 @@ third_party_presolve_status_t third_party_presolve_t::apply_papilo( f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads) + i_t num_cpu_threads, + i_t max_rounds, + i_t max_badgesize) { raft::common::nvtx::range fun_scope("Apply Papilo presolve on host"); @@ -842,11 +855,32 @@ third_party_presolve_status_t third_party_presolve_t::apply_papilo( relative_tolerance, time_limit, dual_postsolve, - num_cpu_threads); - set_presolve_parameters(papilo_presolver, category, original_n_cons, original_n_vars); + num_cpu_threads, + max_rounds); + set_presolve_parameters( + papilo_presolver, category, original_n_cons, original_n_vars, max_badgesize); papilo_presolver.setVerbosityLevel(papilo::VerbosityLevel::kQuiet); - - auto result = papilo_presolver.apply(papilo_problem); + CUOPT_LOG_INFO( + "PRESOLVE_PAPILO_BUDGET rounds=%d badge_cap=%d tlim=%g", max_rounds, max_badgesize, time_limit); + + const auto papilo_t0 = std::chrono::steady_clock::now(); + auto result = papilo_presolver.apply(papilo_problem); + const double papilo_wall = + std::chrono::duration(std::chrono::steady_clock::now() - papilo_t0).count(); + // The effective badge is what set_presolve_parameters actually installed; the cap alone is + // misleading because it only binds once ncols/2 exceeds it. + int effective_badge = std::max(original_n_vars / 2, 32); + if (max_badgesize > 0) { effective_badge = std::min(effective_badge, max_badgesize); } + // hit_tlim distinguishes "presolve converged" from "presolve was cut off mid-round", which + // changes how the reduced problem below should be read. + CUOPT_LOG_INFO( + "PRESOLVE_PAPILO wall=%.3f tlim=%g hit_tlim=%d rounds_cap=%d badge_cap=%d badge_effective=%d", + papilo_wall, + time_limit, + (int)(papilo_wall >= 0.99 * time_limit), + max_rounds, + max_badgesize, + effective_badge); check_presolve_status(result.status); auto status = convert_papilo_presolve_status_to_third_party_presolve_status(result.status); if (result.status == papilo::PresolveStatus::kInfeasible || @@ -902,7 +936,9 @@ third_party_presolve_t::apply_presolve_from_op_problem( f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads) + i_t num_cpu_threads, + i_t max_rounds, + i_t max_badgesize) { auto* handle = op_problem.get_handle_ptr(); @@ -922,7 +958,9 @@ third_party_presolve_t::apply_presolve_from_op_problem( absolute_tolerance, relative_tolerance, time_limit, - num_cpu_threads); + num_cpu_threads, + max_rounds, + max_badgesize); // On terminal statuses the mps entry returns an empty reduced problem; // mirror that shape on the device side without going through H->D. @@ -962,7 +1000,9 @@ third_party_presolve_t::apply_presolve_from_mps_data( f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads) + i_t num_cpu_threads, + i_t max_rounds, + i_t max_badgesize) { presolver_ = presolver; maximize_ = mps.get_sense(); @@ -1012,7 +1052,9 @@ third_party_presolve_t::apply_presolve_from_mps_data( absolute_tolerance, relative_tolerance, time_limit, - num_cpu_threads); + num_cpu_threads, + max_rounds, + max_badgesize); if (status == third_party_presolve_status_t::INFEASIBLE || status == third_party_presolve_status_t::UNBOUNDED || @@ -1074,8 +1116,11 @@ third_party_presolve_status_t third_party_presolve_t::apply_to_subprob settings.dual_tol, time_limit, dual_postsolve, - num_threads); - set_presolve_parameters(papilo_presolver, problem_category_t::MIP, orig_rows, orig_cols); + num_threads, + -1); + // Node presolve already runs under a finite time limit, so it keeps the unbounded round count and + // uncapped badge; the budgets apply to root presolve only. + set_presolve_parameters(papilo_presolver, problem_category_t::MIP, orig_rows, orig_cols, -1); // Disable papilo logs papilo_presolver.setVerbosityLevel(papilo::VerbosityLevel::kQuiet); diff --git a/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp b/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp index 7ed62ef07c..f89812e4eb 100644 --- a/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp +++ b/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp @@ -82,7 +82,9 @@ class third_party_presolve_t { f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads = 0); + i_t num_cpu_threads = 0, + i_t max_rounds = -1, + i_t max_badgesize = -1); // Host entry: takes an mps_data_model_t and returns a host-side reduced // mps_data_model_t. Pure-host throughout @@ -94,7 +96,9 @@ class third_party_presolve_t { f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads = 0); + i_t num_cpu_threads = 0, + i_t max_rounds = -1, + i_t max_badgesize = -1); // Apply the presolve on an simplex::user_problem in-place. Used in sub MIP and (in the future) // restarts. @@ -164,7 +168,9 @@ class third_party_presolve_t { f_t absolute_tolerance, f_t relative_tolerance, double time_limit, - i_t num_cpu_threads); + i_t num_cpu_threads, + i_t max_rounds, + i_t max_badgesize); // Host-only per-backend postsolve helpers. Both resize their vector args // to original-problem dimensions. diff --git a/cpp/src/mip_heuristics/solve.cu b/cpp/src/mip_heuristics/solve.cu index 64d78efbc0..ba041acbd2 100644 --- a/cpp/src/mip_heuristics/solve.cu +++ b/cpp/src/mip_heuristics/solve.cu @@ -13,6 +13,7 @@ #include #include #include +#include #include #include #include @@ -332,6 +333,36 @@ mip_solution_t run_mip_solver( } } +namespace mip { + +// Features of the problem as the user handed it in, i.e. before any reduction. This is what Papilo +// itself will work on, so its budget is derived from these rather than from the reduced problem. +template +presolve_features_t papilo_presolve_features(optimization_problem_t const& op_problem) +{ + presolve_features_t f{}; + f.n_vars = op_problem.get_n_variables(); + f.n_cons = op_problem.get_n_constraints(); + f.nnz = op_problem.get_nnz(); + + const auto var_types = op_problem.get_variable_types_host(); + const auto lower = op_problem.get_variable_lower_bounds_host(); + const auto upper = op_problem.get_variable_upper_bounds_host(); + for (size_t j = 0; j < var_types.size(); ++j) { + if (var_types[j] != var_t::INTEGER) { continue; } + f.n_int += 1.0; + if (lower[j] >= 0.0 && upper[j] <= 1.0) { f.n_bin += 1.0; } + } + + const auto offsets = op_problem.get_constraint_matrix_offsets_host(); + for (size_t i = 0; i + 1 < offsets.size(); ++i) { + f.max_row_len = std::max(f.max_row_len, offsets[i + 1] - offsets[i]); + } + return f; +} + +} // namespace mip + template mip_solution_t solve_mip_helper(optimization_problem_t& op_problem, mip_solver_settings_t const& settings_const) @@ -561,12 +592,19 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p sort_csr(op_problem); // allocate not more than 10% of the time limit to presolve. // Note that this is not the presolve time, but the time limit for presolve. - const auto& hp = settings.heuristic_params; - double presolve_time_limit = - std::min(hp.presolve_time_ratio * time_limit, hp.presolve_max_time); + const auto& hp = settings.heuristic_params; + const auto papilo_features = mip::papilo_presolve_features(op_problem); + const auto papilo_budget = mip::evaluate_presolve_budget(hp, papilo_features); + mip::log_presolve_budget("PAPILO", hp.presolve_budget_policy, papilo_features, papilo_budget); + + // The round and badge caps are what shape presolve; the timer only carries the time actually + // left in the solve, so it cannot overrun the whole budget but also does not truncate + // presolve at some fraction of it the way a presolve-specific cap did. + double presolve_time_limit = timer.remaining_time(); if (settings.determinism_mode == CUOPT_MODE_DETERMINISTIC) { presolve_time_limit = std::numeric_limits::infinity(); } + presolver = std::make_unique>(); auto result = presolver->apply_presolve_from_op_problem( op_problem, @@ -576,7 +614,9 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p settings.tolerances.absolute_tolerance, settings.tolerances.relative_tolerance, presolve_time_limit, - settings.num_cpu_threads); + settings.num_cpu_threads, + papilo_budget.papilo_max_rounds, + papilo_budget.papilo_max_badgesize); if (result.status == mip::third_party_presolve_status_t::INFEASIBLE) { return mip_solution_t(mip_termination_status_t::Infeasible, @@ -606,6 +646,19 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p CUOPT_LOG_INFO("%d implied integers", presolve_result_opt->implied_integer_indices.size()); } CUOPT_LOG_INFO("Papilo presolve time: %.2f", presolve_time); + // What the round cap actually bought, logged here rather than inferred from the probing stage + // so it is still recorded when the run never gets that far. + CUOPT_LOG_INFO( + "PRESOLVE_PAPILO_REDUCED nvars=%d ncons=%d nnz=%d nint=%d nbin=%d from_nvars=%.0f " + "from_ncons=%.0f from_nnz=%.0f", + problem.n_variables, + problem.n_constraints, + problem.nnz, + problem.n_integer_vars, + problem.n_binary_vars, + papilo_features.n_vars, + papilo_features.n_cons, + papilo_features.nnz); if (result.status == mip::third_party_presolve_status_t::OPTIMAL) { CUOPT_LOG_INFO("Optimal solution found during presolve."); From a74d4860efdb190e5b5b2d29a317de256bfebd9b Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 11:15:27 +0200 Subject: [PATCH 02/20] fix(mip): scale the probing budget with the candidate count The probing-cache budget was an absolute number of work units, and every derived policy was additionally clamped to at most 120. Probing cost is close to linear in the candidate count, so a constant budget probes a small problem exhaustively and a large one barely at all: at the default of 30 units that was 21% of 30n20b8 but 0.8% of netdiversion. On the 240-instance benchmark this cost four instances their feasible solution (physiciansched3-3, neos-3216931-puriri, 30n20b8, neos-5104907-jarama), which is most of the 10.97 -> 12.33 mean error regression. 169 of 239 instances exhausted the budget, mean error delta +2.18, against -0.59 for the 70 that did not. None of the 169 were limited by the clock: probing stopped after 0.2-1.8s with ~590s still available. Policies now state the fraction of candidates they want probed, converted with the same cost model the probing loop charges. A fraction >= 1 is carried as infinity rather than a large finite number, since the per-candidate cost is an average and a budget sized for full coverage would still truncate an instance whose probes are dearer than average. probing_step_size drops from 512 to 128, because at 512 the first step already overshot a 30-unit budget by 2-2.7x and probe counts were pinned to multiples of the step. `fixed` becomes the Papilo-only arm (probing unbounded, rounds and badge capped), and `size` keeps the measured badge rule at a quarter coverage, so benchmarking the two separates the Papilo caps from the probing budget. Analysis and numbers in design_summaries/presolve_budget/. Signed-off-by: akif --- .../presolve/presolve_budget_policy.hpp | 68 +++++++++++++------ 1 file changed, 47 insertions(+), 21 deletions(-) diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index 891b598cb2..c09ed75da1 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -76,6 +76,9 @@ struct presolve_budget_t { double probing_work_limit{std::numeric_limits::infinity()}; // Probed variables per step, i.e. the granularity at which the budget can be enforced. int probing_step_size{2048}; + // Coverage the policy asked for, kept only so the log can be compared against what was realised; + // the two diverge when an instance's probes are dearer than the average the budget assumed. + double intended_probe_fraction{1.0}; }; namespace detail { @@ -98,15 +101,25 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(feat.nnz, 1.0); - const double arl = std::max(feat.avg_row_len(), 1.0); - const double n_int = std::max(feat.n_int, 1.0); - const double n_bin = std::max(feat.n_bin, 1.0); - const double bf = feat.bin_frac(); - + const double nnz = std::max(feat.nnz, 1.0); + const double arl = std::max(feat.avg_row_len(), 1.0); + // Probing candidates are the integers of the problem the probing cache runs on. + const double n_cand = std::max(feat.n_int, 1.0); + const double n_bin = std::max(feat.n_bin, 1.0); + const double bf = feat.bin_frac(); + + // Probing cost is close to linear in the candidate count -- ~0.055 work units per candidate, + // measured across the 240-instance benchmark -- and the candidate set is the integers of the + // reduced problem. A constant budget therefore probes a small problem exhaustively and a large + // one barely at all: a flat 30 units covered 21% of 30n20b8 but 0.8% of netdiversion, which cost + // four instances their feasible solution. Policies state the fraction of candidates they want + // probed instead, converted below with the same cost model the probing loop charges. + constexpr double avg_iters_per_probe = 3.5; + const double per_candidate_work = + (double)hp.probe_host_overhead_work + avg_iters_per_probe * (double)hp.probe_iter_work; // Lets a single knob scale every derived policy without recompiling: 1.0 at the default of 30. const double work_scale = static_cast(hp.cuopt_presolve_work_limit) / 30.0; - double raw_work = 30.0; + double probe_fraction = 1.0; switch (static_cast(hp.presolve_budget_policy)) { case presolve_budget_policy_t::legacy: @@ -116,12 +129,14 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t 2.0e4; b.papilo_max_rounds = wide ? 50 : 20; b.papilo_max_badgesize = wide ? 32 : -1; - raw_work = 120.0; - b.probing_step_size = 512; + probe_fraction = 0.25; + b.probing_step_size = 128; break; } @@ -151,8 +170,8 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t= 0.9 ? 50 : bf >= 0.5 ? 30 : 15; b.papilo_max_badgesize = clamp_i(std::max(n_bin / 2.0, 32.0), 32, 1024); - raw_work = (2.0 / 3.0) * std::sqrt(n_bin); - b.probing_step_size = 512; + probe_fraction = clamp_d(bf, 0.05, 1.0); + b.probing_step_size = 128; break; // Multiplicative over the three effects above, so no single feature can dominate the budget. @@ -172,13 +191,19 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t= 1.0 ? std::numeric_limits::infinity() + : probe_fraction * n_cand * per_candidate_work; return b; } @@ -192,7 +217,7 @@ inline void log_presolve_budget(const char* stage, CUOPT_LOG_INFO( "PRESOLVE_BUDGET stage=%s policy=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f nbin=%.0f " "arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " - "rounds=%d badge=%d work=%.3f step=%d", + "rounds=%d badge=%d work=%.3f step=%d intended_probe_frac=%.4f", stage, presolve_budget_policy_name(policy), f.n_vars, @@ -209,7 +234,8 @@ inline void log_presolve_budget(const char* stage, b.papilo_max_rounds, b.papilo_max_badgesize, b.probing_work_limit, - b.probing_step_size); + b.probing_step_size, + b.intended_probe_fraction); } } // namespace cuopt::mathematical_optimization::mip From ca84fecb3fbcc9534d970bca81d95c9b8ae2f882 Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 11:31:41 +0200 Subject: [PATCH 03/20] feat(mip): select presolve budget points with CUOPT_CONFIG_ID The 240-instance run could not separate the Papilo round/badge caps from the probing coverage, because a single policy fixes both at once. Add a 2x3 factorial over the two -- two Papilo rules against three probing coverages -- selected by CUOPT_CONFIG_ID so one build covers the whole sweep: 0/1/2 fixed (rounds 30, badge 1024) x 1.00 / 0.25 / 0.05 3/4/5 size (wide 50/32, narrow 20/uncapped) x 1.00 / 0.25 / 0.05 Configs 0 and 3 leave probing unbounded and so isolate the Papilo rule. Coverage is spaced geometrically rather than evenly because truncating probing won at 0.5-4.3% of candidates (bab6, square41, square47) and lost at 9.1-21.6% (physiciansched3-3, 30n20b8), which puts the crossover below 10%. The id is resolved once inside evaluate_presolve_budget rather than at each call site, so the Papilo stage and the probing stage cannot end up running different points of the sweep. It overrides the policy hyper-parameter, and the effective policy and config now travel on presolve_budget_t so PRESOLVE_BUDGET stays attributable instead of reporting the hyper-parameter that was overridden. Verified on 30n20b8: realised coverage 100 / 24.3 / 5.4% for configs 0-2 and 100 / 24.9 / 5.5% for 3-5. Signed-off-by: akif --- .../diversity/diversity_manager.cu | 2 +- .../presolve/presolve_budget_policy.hpp | 72 +++++++++++++++++-- cpp/src/mip_heuristics/solve.cu | 2 +- 3 files changed, 69 insertions(+), 7 deletions(-) diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 02b3fc4c62..5803c20440 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -317,7 +317,7 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ run_probing_cache = false; } if (run_probing_cache) { - log_presolve_budget("PROBING", hp.presolve_budget_policy, probing_features, probing_budget); + log_presolve_budget("PROBING", probing_features, probing_budget); // Run probing cache before trivial presolve to discover variable implications. The work budget // decides how much probing happens; the timer only carries the time actually left in the solve, // so it cannot overrun the whole budget but also is not what shapes the amount of probing. diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index c09ed75da1..dd8f56b458 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -13,7 +13,10 @@ #include #include +#include +#include #include +#include namespace cuopt::mathematical_optimization::mip { @@ -46,6 +49,51 @@ inline const char* presolve_budget_policy_name(int policy) } } +// Benchmark points selected by CUOPT_CONFIG_ID, overriding the policy hyper-parameter so one build +// covers the whole sweep. The 240-instance run conflated two effects -- the Papilo round/badge caps +// and how much of the problem gets probed -- so these span them as a 2x3 factorial: two Papilo +// rules against three probing coverages. Coverage is spaced geometrically rather than evenly +// because truncating probing both won (bab6 at 0.5% of candidates, square41 at 4.3%) and lost +// (30n20b8 at 21.6%, physiciansched3-3 at 9.1%), which puts the crossover somewhere below 10%. +struct presolve_config_t { + presolve_budget_policy_t papilo_rule; + double probe_fraction; +}; + +inline constexpr presolve_config_t presolve_configs[] = { + {presolve_budget_policy_t::fixed, 1.00}, // 0: Papilo rounds=30/badge=1024, probing unbounded + {presolve_budget_policy_t::fixed, 0.25}, // 1 + {presolve_budget_policy_t::fixed, 0.05}, // 2 + {presolve_budget_policy_t::size, 1.00}, // 3: Papilo measured wide/narrow rule, probing unbounded + {presolve_budget_policy_t::size, 0.25}, // 4 + {presolve_budget_policy_t::size, 0.05}, // 5 +}; +inline constexpr int n_presolve_configs = 6; + +// -1 when unset or out of range. Read once; the environment cannot change mid-run. +inline int presolve_config_id() +{ + static const int id = []() -> int { + const char* raw = std::getenv("CUOPT_CONFIG_ID"); + if (raw == nullptr) { return -1; } + try { + const int v = std::stoi(raw); + if (v < 0 || v >= n_presolve_configs) { + CUOPT_LOG_WARN("CUOPT_CONFIG_ID=%d is outside [0, %d); ignoring it for presolve budgets", + v, + n_presolve_configs); + return -1; + } + CUOPT_LOG_INFO("Using presolve budget config %d from CUOPT_CONFIG_ID", v); + return v; + } catch (const std::exception& e) { + CUOPT_LOG_WARN("Failed to parse CUOPT_CONFIG_ID: %s", e.what()); + return -1; + } + }(); + return id; +} + // Dimensions plus cheap structural ratios. Both presolve stages populate this from whatever problem // representation they hold: Papilo from the original problem before any reduction, the cuOpt // probing cache from the Papilo-reduced problem. The two therefore see different feature values for @@ -79,6 +127,10 @@ struct presolve_budget_t { // Coverage the policy asked for, kept only so the log can be compared against what was realised; // the two diverge when an instance's probes are dearer than the average the budget assumed. double intended_probe_fraction{1.0}; + // The policy that actually ran and the config that selected it, so a log line is attributable + // even when CUOPT_CONFIG_ID overrode the hyper-parameter. + int policy{1}; + int config_id{-1}; }; namespace detail { @@ -121,7 +173,16 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(hp.cuopt_presolve_work_limit) / 30.0; double probe_fraction = 1.0; - switch (static_cast(hp.presolve_budget_policy)) { + // A config id, when set, replaces both the policy and its coverage. Resolved here rather than at + // each call site so the Papilo stage and the probing stage cannot disagree about which point of + // the sweep is running. + const int config = presolve_config_id(); + const auto policy = config >= 0 ? presolve_configs[config].papilo_rule + : (presolve_budget_policy_t)hp.presolve_budget_policy; + b.policy = (int)policy; + b.config_id = config; + + switch (policy) { case presolve_budget_policy_t::legacy: b.papilo_max_rounds = -1; b.papilo_max_badgesize = -1; @@ -200,6 +261,7 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t= 0) { probe_fraction = presolve_configs[config].probe_fraction; } probe_fraction *= work_scale; b.intended_probe_fraction = std::min(probe_fraction, 1.0); b.probing_work_limit = probe_fraction >= 1.0 ? std::numeric_limits::infinity() @@ -210,16 +272,16 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t solve_mip_helper(optimization_problem_t& op_p const auto& hp = settings.heuristic_params; const auto papilo_features = mip::papilo_presolve_features(op_problem); const auto papilo_budget = mip::evaluate_presolve_budget(hp, papilo_features); - mip::log_presolve_budget("PAPILO", hp.presolve_budget_policy, papilo_features, papilo_budget); + mip::log_presolve_budget("PAPILO", papilo_features, papilo_budget); // The round and badge caps are what shape presolve; the timer only carries the time actually // left in the solve, so it cannot overrun the whole budget but also does not truncate From 0360daac723b6047c2f3922a288aaa05d55eddf1 Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 13:35:00 +0200 Subject: [PATCH 04/20] fix(mip): keep a wall ceiling on Papilo and probing Commit 8de51122f made both presolve timers carry remaining_time(), dropping min(presolve_time_ratio * time_limit, presolve_max_time) for Papilo and max_time_on_probing for the probing cache. Those were the only bounds on presolve in time, and a work budget does not substitute for one: realised probing throughput spans 12 to 450 work units per second, so on s100 even a 25% coverage budget ran 583-593s of a 600s solve. ns1760995 converges in well under 30 rounds, so its round cap never binds and Papilo ran 243-403s in every config. Presolve starving the solve is what drove the inf mip_gaps up: 16 of the 21 across the six configs have presolve taking 40-99.8% of the budget, and the probing-unbounded configs show 7 each against 1-3 for the truncated ones. With no time left the root LP never produces a dual bound. Restore both as min(cap, remaining_time()), so the cap is a hard ceiling and remaining_time only stops presolve reaching past the end of the solve. The structural budgets still shape presolve; the ceiling only stops it starving the solve. Verified at a 180s limit: ns1760995 Papilo 18.2s against tlim=18 with hit_tlim=1, supportcase6 probing 60.2s against wall_limit=60. Signed-off-by: akif --- .../mip_heuristics/diversity/diversity_manager.cu | 11 ++++++++--- cpp/src/mip_heuristics/solve.cu | 12 ++++++++---- 2 files changed, 16 insertions(+), 7 deletions(-) diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 5803c20440..a7e98e8773 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -319,9 +319,14 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ if (run_probing_cache) { log_presolve_budget("PROBING", probing_features, probing_budget); // Run probing cache before trivial presolve to discover variable implications. The work budget - // decides how much probing happens; the timer only carries the time actually left in the solve, - // so it cannot overrun the whole budget but also is not what shapes the amount of probing. - timer_t probing_timer{global_timer.remaining_time()}; + // shapes how much probing happens but cannot bound how long that takes: realised throughput + // spans 12 to 450 work units per second across the benchmark, so at the slow end even a quarter + // of the candidates ran 593 of the 600 available seconds and left the root LP with no dual + // bound. max_time_on_probing is therefore a hard ceiling rather than a target, and + // remaining_time only stops it reaching past the end of the solve. + f_t time_for_probing_cache = + std::min((f_t)diversity_config.max_time_on_probing, (f_t)global_timer.remaining_time()); + timer_t probing_timer{time_for_probing_cache}; const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem bool problem_is_infeasible = compute_probing_cache(ls.constraint_prop.bounds_update, diff --git a/cpp/src/mip_heuristics/solve.cu b/cpp/src/mip_heuristics/solve.cu index a5a3950dfb..1511ce1c0e 100644 --- a/cpp/src/mip_heuristics/solve.cu +++ b/cpp/src/mip_heuristics/solve.cu @@ -597,10 +597,14 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p const auto papilo_budget = mip::evaluate_presolve_budget(hp, papilo_features); mip::log_presolve_budget("PAPILO", papilo_features, papilo_budget); - // The round and badge caps are what shape presolve; the timer only carries the time actually - // left in the solve, so it cannot overrun the whole budget but also does not truncate - // presolve at some fraction of it the way a presolve-specific cap did. - double presolve_time_limit = timer.remaining_time(); + // The round and badge caps shape presolve, but they do not bound its cost: ns1760995 + // converges in well under 30 rounds and still spent 243-403s of a 600s budget in Papilo, + // which left no time to find a dual bound. The wall cap is what makes that survivable, so it + // stays a hard ceiling and remaining_time only stops presolve reaching past the end of the + // solve. + double presolve_time_limit = + std::min(std::min(hp.presolve_time_ratio * time_limit, hp.presolve_max_time), + timer.remaining_time()); if (settings.determinism_mode == CUOPT_MODE_DETERMINISTIC) { presolve_time_limit = std::numeric_limits::infinity(); } From 568da86b18e21529a7d3b449f38f9c54bac2bac5 Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 17:35:49 +0200 Subject: [PATCH 05/20] feat(mip): bound presolve by cost proxies so it concludes before its ceiling Restoring the wall ceilings stopped presolve starving the solve, but only by truncating mid-round, which pays for reduction and then discards it. Size each stage's own limit so it finishes well inside 60s and the ceiling never fires. Probing throughput spans 1.9 to 689 work units/s over 660 measured runs, so no coverage fraction can bound time. nnz + n_cand * avg_col_len tracks the cost of a propagation sweep best of the proxies tried, and a ceiling of 1.5e8 over it caps the slowest measured run at 44s with none above 60s, turning today's 597s runaways into 3-34s. For Papilo, rounds never bind: all 240 instances converged inside 30 without hitting a round or time cap, and minbadgesize is the whole cost. Clamp it to 32 above n_bin * avg_col_len of 2e5, where a large badge stops paying for itself -- square47 40.6s to 8.1s and sorrell3 42.9s to 14.0s, both bit-identical reduced problems. ns1760995 is the exception where the badge does buy reduction, but only at 1024 and at 274s, and that is what loses the instance: truncating Papilo at 60s closed it to 1.37% while every run that let it finish ended with no dual bound. --- .../mathematical_optimization/constants.h | 4 +- .../mip/heuristics_hyper_params.hpp | 8 +- cpp/src/math_optimization/solver_settings.cu | 3 +- .../presolve/presolve_budget_policy.hpp | 93 +++++++++++++++---- 4 files changed, 87 insertions(+), 21 deletions(-) diff --git a/cpp/include/cuopt/mathematical_optimization/constants.h b/cpp/include/cuopt/mathematical_optimization/constants.h index da92ac8733..6ef2de9cee 100644 --- a/cpp/include/cuopt/mathematical_optimization/constants.h +++ b/cpp/include/cuopt/mathematical_optimization/constants.h @@ -109,7 +109,9 @@ #define CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE "mip_hyper_heuristic_probing_step_size" #define CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK \ "mip_hyper_heuristic_probe_host_overhead_work" -#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK "mip_hyper_heuristic_probe_iter_work" +#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK "mip_hyper_heuristic_probe_iter_work" +#define CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE \ + "mip_hyper_heuristic_probing_work_time_scale" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO "mip_hyper_heuristic_root_lp_time_ratio" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME "mip_hyper_heuristic_root_lp_max_time" #define CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT "mip_hyper_heuristic_rins_time_limit" diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index e2ed78e2e4..5c0d6ebbc8 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -27,7 +27,7 @@ struct mip_heuristics_hyper_params_t { // Presolve budgeting. presolve_budget_policy selects how the four knobs below are derived from // the problem's dimensions and structure (see presolve_budget_policy.hpp); the values here are // the defaults the policy starts from and the literal values used by the `manual` policy. - i_t presolve_budget_policy = 1; // presolve_budget_policy_t + i_t presolve_budget_policy = 7; // presolve_budget_policy_t (cost) i_t presolve_max_rounds = 30; // Papilo presolve rounds cap (<=0 = Papilo default) i_t papilo_probing_max_badgesize = 1024; // ceiling on Papilo's probing.minbadgesize f_t cuopt_presolve_work_limit = 30.0; // probing-cache budget, work units @@ -37,6 +37,12 @@ struct mip_heuristics_hyper_params_t { // differs between instances. f_t probe_host_overhead_work = 0.02; // charged per probed variable f_t probe_iter_work = 0.01; // charged per multi-probe propagation iteration + // Numerator of the work ceiling that keeps probing inside a wall-clock budget. Work units are + // reproducible but their cost is not: realised throughput spans 1.9 to 689 units/s across the + // 240-instance benchmark, so a coverage target alone lets slow instances run for the whole solve. + // Dividing this by the cost proxy (nnz + n_cand * avg_col_len) bounds the wall time instead; at + // 1.5e8 the slowest of 660 measured probing runs finishes in 44s, none above 60s. + f_t probing_work_time_scale = 1.5e8; f_t root_lp_time_ratio = 0.1; // fraction of total time for root LP f_t root_lp_max_time = 15.0; // hard cap on root LP seconds diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index bedc60a053..370cd0b3f6 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -109,6 +109,7 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings {CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT, &mip_settings.heuristic_params.cuopt_presolve_work_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(30.0), "probing-cache budget in work units (manual policy, or ceiling for derived policies)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK, &mip_settings.heuristic_params.probe_host_overhead_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.02), "work units charged per probed variable (host overhead)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK, &mip_settings.heuristic_params.probe_iter_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.01), "work units charged per multi-probe propagation iteration"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(1.5e8), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, @@ -170,7 +171,7 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") {CUOPT_MIP_HYPER_HEURISTIC_POPULATION_SIZE, &mip_settings.heuristic_params.population_size, 1, std::numeric_limits::max(), 32, "max solutions in pool"}, {CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS, &mip_settings.heuristic_params.num_cpufj_threads, 0, std::numeric_limits::max(), 8, "parallel CPU FJ climbers"}, - {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY, &mip_settings.heuristic_params.presolve_budget_policy, 0, 6, 1, "how presolve budgets are derived: 0=legacy 1=fixed 2=size 3=density 4=binary 5=combined 6=manual"}, + {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY, &mip_settings.heuristic_params.presolve_budget_policy, 0, 7, 7, "how presolve budgets are derived: 0=legacy 1=fixed 2=size 3=density 4=binary 5=combined 6=manual 7=cost"}, {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS, &mip_settings.heuristic_params.presolve_max_rounds, -1, std::numeric_limits::max(), 30, "Papilo presolve rounds cap (<=0 keeps Papilo default)"}, {CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE, &mip_settings.heuristic_params.papilo_probing_max_badgesize, -1, std::numeric_limits::max(), 1024, "ceiling on Papilo probing.minbadgesize (<=0 leaves it uncapped)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE, &mip_settings.heuristic_params.probing_step_size, 1, std::numeric_limits::max(), 512, "probed variables between probing-cache work-budget checks"}, diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index dd8f56b458..724a6c37ef 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -33,6 +33,7 @@ enum class presolve_budget_policy_t : int { binary = 4, combined = 5, manual = 6, + cost = 7, }; inline const char* presolve_budget_policy_name(int policy) @@ -45,28 +46,39 @@ inline const char* presolve_budget_policy_name(int policy) case presolve_budget_policy_t::binary: return "binary"; case presolve_budget_policy_t::combined: return "combined"; case presolve_budget_policy_t::manual: return "manual"; + case presolve_budget_policy_t::cost: return "cost"; default: return "unknown"; } } // Benchmark points selected by CUOPT_CONFIG_ID, overriding the policy hyper-parameter so one build -// covers the whole sweep. The 240-instance run conflated two effects -- the Papilo round/badge caps -// and how much of the problem gets probed -- so these span them as a 2x3 factorial: two Papilo -// rules against three probing coverages. Coverage is spaced geometrically rather than evenly -// because truncating probing both won (bab6 at 0.5% of candidates, square41 at 4.3%) and lost -// (30n20b8 at 21.6%, physiciansched3-3 at 9.1%), which puts the crossover somewhere below 10%. +// covers the whole sweep. The previous 240-instance run showed that neither stage was actually +// bounded: presolve ran 243-403s (Papilo, ns1760995) and up to 598s (probing, six instances) out of +// a 600s solve, and the instances that lost their dual bound are the ones that spent it there. Both +// stages now carry a cost model instead of a coverage target, so these points span the two models +// against the previously measured coverage fractions. +// +// `work_time_scale` of 0 disables the probing ceiling, which separates "probe a fixed fraction" +// from "probe until a wall-clock proxy is exhausted"; the fractions bracket the crossover the +// earlier run pointed at, where truncation won on bab6 (0.5% of candidates), square41 (4.3%) and +// square47 (2.9%) and lost on 30n20b8 (21.6%) and physiciansched3-3 (9.1%). struct presolve_config_t { presolve_budget_policy_t papilo_rule; double probe_fraction; + double work_time_scale; }; inline constexpr presolve_config_t presolve_configs[] = { - {presolve_budget_policy_t::fixed, 1.00}, // 0: Papilo rounds=30/badge=1024, probing unbounded - {presolve_budget_policy_t::fixed, 0.25}, // 1 - {presolve_budget_policy_t::fixed, 0.05}, // 2 - {presolve_budget_policy_t::size, 1.00}, // 3: Papilo measured wide/narrow rule, probing unbounded - {presolve_budget_policy_t::size, 0.25}, // 4 - {presolve_budget_policy_t::size, 0.05}, // 5 + // Old Papilo rule, so the probing ceiling is measured on its own. + {presolve_budget_policy_t::fixed, 1.00, 1.5e8}, // 0 + // New Papilo cost rule against the ceiling and the two coverage fractions. + {presolve_budget_policy_t::cost, 1.00, 1.5e8}, // 1 + {presolve_budget_policy_t::cost, 0.25, 1.5e8}, // 2 + {presolve_budget_policy_t::cost, 0.05, 1.5e8}, // 3 + // Half the ceiling: 29s worst case rather than 44s, to see whether the margin costs quality. + {presolve_budget_policy_t::cost, 1.00, 7.5e7}, // 4 + // Coverage only, no ceiling -- the control that says whether the ceiling earns its keep. + {presolve_budget_policy_t::cost, 0.05, 0.0}, // 5 }; inline constexpr int n_presolve_configs = 6; @@ -127,6 +139,10 @@ struct presolve_budget_t { // Coverage the policy asked for, kept only so the log can be compared against what was realised; // the two diverge when an instance's probes are dearer than the average the budget assumed. double intended_probe_fraction{1.0}; + // The work ceiling the cost proxy produced, and whether it rather than the coverage target is + // what set probing_work_limit. Logged so a run can be attributed to one or the other offline. + double probing_work_ceiling{std::numeric_limits::infinity()}; + bool probing_ceiling_binding{false}; // The policy that actually ran and the config that selected it, so a log line is attributable // even when CUOPT_CONFIG_ID overrode the hyper-parameter. int policy{1}; @@ -172,6 +188,14 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(hp.cuopt_presolve_work_limit) / 30.0; double probe_fraction = 1.0; + double work_time_scale = (double)hp.probing_work_time_scale; + + // Cost of one probing sweep: every propagation touches the rows of the probed column, so the work + // a second buys falls off with problem size. nnz + n_cand * avg_col_len tracked that better than + // nnz, arl or n_cand alone over 660 measured probing runs, and dividing the scale by it bounds + // the wall time that a pure coverage target cannot: throughput ranged 1.9 to 689 work units per + // second, so the same budget was worth 360x more time on one instance than another. + const double probing_cost_proxy = nnz + n_cand * std::max(feat.avg_col_len(), 1.0); // A config id, when set, replaces both the policy and its coverage. Resolved here rather than at // each call site so the Papilo stage and the probing stage cannot disagree about which point of @@ -225,6 +249,26 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(feat.avg_col_len(), 1.0); + b.papilo_max_rounds = 30; + b.papilo_max_badgesize = papilo_probe_cost > 2.0e5 ? 32 : -1; + probe_fraction = 1.0; + b.probing_step_size = 128; + break; + } + // A single propagation sweep costs roughly one pass over each touched row, so long rows make // every probe more expensive: buy fewer of them, and size the badge so one badge's working // limit (~2*nnz in Papilo) stays roughly constant. @@ -258,14 +302,24 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t= 0) { + probe_fraction = presolve_configs[config].probe_fraction; + work_time_scale = presolve_configs[config].work_time_scale; + } + probe_fraction *= work_scale; + b.intended_probe_fraction = std::min(probe_fraction, 1.0); + // A fraction at or above 1 means "probe everything", carried as no budget rather than as a large // number: per_candidate_work is an average, so on an instance whose probes are dearer than // average a finite budget sized for full coverage would still cut probing short. - if (config >= 0) { probe_fraction = presolve_configs[config].probe_fraction; } - probe_fraction *= work_scale; - b.intended_probe_fraction = std::min(probe_fraction, 1.0); - b.probing_work_limit = probe_fraction >= 1.0 ? std::numeric_limits::infinity() - : probe_fraction * n_cand * per_candidate_work; + const double coverage_work = probe_fraction >= 1.0 ? std::numeric_limits::infinity() + : probe_fraction * n_cand * per_candidate_work; + b.probing_work_ceiling = work_time_scale > 0.0 ? work_time_scale / probing_cost_proxy + : std::numeric_limits::infinity(); + // Whichever binds first. The ceiling is what keeps probing away from the wall cap it used to run + // into, and the fraction is what stops it probing more of a cheap problem than is worth probing. + b.probing_ceiling_binding = b.probing_work_ceiling < coverage_work; + b.probing_work_limit = std::min(coverage_work, b.probing_work_ceiling); return b; } @@ -278,7 +332,8 @@ inline void log_presolve_budget(const char* stage, CUOPT_LOG_INFO( "PRESOLVE_BUDGET stage=%s config=%d policy=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f " "nbin=%.0f arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " - "rounds=%d badge=%d work=%.3f step=%d intended_probe_frac=%.4f", + "rounds=%d badge=%d work=%.3f step=%d intended_probe_frac=%.4f work_ceiling=%.3f " + "ceiling_binding=%d", stage, b.config_id, presolve_budget_policy_name(b.policy), @@ -297,7 +352,9 @@ inline void log_presolve_budget(const char* stage, b.papilo_max_badgesize, b.probing_work_limit, b.probing_step_size, - b.intended_probe_fraction); + b.intended_probe_fraction, + b.probing_work_ceiling, + (int)b.probing_ceiling_binding); } } // namespace cuopt::mathematical_optimization::mip From e3226169784c17ab93280851324ee3d1f87604a4 Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 17:44:37 +0200 Subject: [PATCH 06/20] cuopt-skill-evolution: add stage-budget and regression-attribution methodology to cuopt-developer Bounding MIP presolve by work units regressed the benchmark, and three separate assumptions had to be falsified by measurement before the fix landed: that a deterministic work counter bounds wall time, that the plausible-looking knob (Papilo rounds) drives cost, and that reduced probing coverage caused a given instance's feasibility loss. The generalizable part is not the presolve numbers but the method: measure the throughput spread before trusting a work budget, size a limit so the stage concludes before its ceiling rather than being truncated by it, check whether a cap ever binds before tuning it, and run the ordered attribution checks -- did the limit bind, did the stage starve a later one, does the failure reproduce unbounded -- before changing a limit's value. --- skills/cuopt-developer/SKILL.md | 8 ++ .../references/stage_budgets.md | 78 +++++++++++++++++++ 2 files changed, 86 insertions(+) create mode 100644 skills/cuopt-developer/references/stage_budgets.md diff --git a/skills/cuopt-developer/SKILL.md b/skills/cuopt-developer/SKILL.md index 5e623e7836..59501a8477 100644 --- a/skills/cuopt-developer/SKILL.md +++ b/skills/cuopt-developer/SKILL.md @@ -265,6 +265,14 @@ When implementing or debugging **VRP dimensions** (constraints, objectives, forw Read it **before** adding a new dimension or changing combine semantics. +## Budgeting a solver stage and attributing regressions + +When adding or tuning a **limit on a stage that shares the solve's time budget** (presolve, probing, cut generation, heuristics), or when a benchmark regression looks like it came from one, read: + +- **`references/stage_budgets.md`** — why a deterministic work counter does not bound wall time, how to size a limit so the stage concludes before its ceiling, and the ordered checks that attribute a regression to a budget before tuning its value. + +Run the attribution checks **before** changing a limit's value — the common failure is tuning a limit that was never what bound. + ## Numerical issues in non-routing solver internals When a bug surfaces as **wrong-but-plausible** solver output (invalid lower bound, unexpectedly large duals, 10× iteration blow-up after a small change) rather than a crash, read: diff --git a/skills/cuopt-developer/references/stage_budgets.md b/skills/cuopt-developer/references/stage_budgets.md new file mode 100644 index 0000000000..7858d81496 --- /dev/null +++ b/skills/cuopt-developer/references/stage_budgets.md @@ -0,0 +1,78 @@ +# Budgeting a solver stage, and attributing a benchmark regression to one + +Applies to any bounded stage that consumes a shared time budget — presolve, probing, cut +generation, root LP, heuristics. + +## A work counter buys reproducibility, not a time bound + +A deterministic work counter (units of effort charged by a cost model) makes a stage's +*coverage* machine-independent. It does not bound how long that coverage takes, because the +work-per-second a stage achieves varies with problem structure. Measure the spread before +assuming otherwise: log realised `work / wall` per run and look at the range. When it spans +orders of magnitude, a single work limit is worth wildly different amounts of time per +instance, and the slow end will consume the whole solve. + +Two limits, two jobs — keep both: + +- **Work limit** — makes coverage reproducible across machines. +- **Wall ceiling** — keeps the stage from starving whatever runs after it. + +## Size the limit to conclude *before* the ceiling + +A stage that hits its wall ceiling is the worst case, not the safe one: the effort is spent +and the partial result is discarded. Prefer a limit the stage finishes inside, leaving the +ceiling as a backstop that never fires. Verify by checking that the stage's `hit_limit` / +`timed_out` flag is clear on the instances that used to be truncated. + +To bound wall time structurally, divide a scale constant by a **cost proxy** built from +dimensions the stage's inner loop actually touches. Rank candidate proxies by how nearly +`throughput * proxy` is constant across measured runs, and pick the scale so the slowest +measured run still finishes with margin. The proxy will be too noisy to *predict* time; it +only needs to *bound* it, so tune it against the worst case and accept that typical instances +get less than they could use. + +## Verify which knob costs, rather than the one that looks like it should + +Iteration or round counts read like cost limits and often are not: a stage that converges +before its round cap is not bounded by it at all. Check whether the cap ever binds +(`rounds_used < rounds_cap`, `hit_tlim=0`) before tuning it. The knob that drives cost is +usually the one controlling work *per* round — batch or candidate-set size — and it may sit +behind a third-party default. + +A knob can also be free on some instances and a genuine quality lever on others. Sweep it and +record both cost and the stage's own output (reduced dimensions, cuts kept) — identical output +at lower cost means the knob is free; changed output means there is a real trade-off, which +must then be settled end-to-end rather than on the stage's own metric. + +## More effort in a stage can be worse end-to-end + +A stage's own metric improving is not evidence the solve improves. Effort spent there is taken +from everything downstream, so compare final objective and gap, not reduction counts. Watch +for a stage consuming a large share of the budget alongside a missing dual bound (infinite +gap): that pairing means the stage starved the root relaxation. + +## Attributing a regression before tuning + +Run these checks in order; each can end the investigation: + +1. **Did the limit bind?** Split instances by the `budget_exhausted` flag and compare the mean + error delta per group. If the regressed set is not the bound set, the limit is not the + cause. +2. **Did the stage starve a later one?** Correlate stage wall time as a fraction of the budget + against instances that lost their dual bound. +3. **Control: remove the limit on the failing instance.** If the failure reproduces with the + stage fully unbounded, the limit is exonerated — look at whichever stage the time actually + went to. This is the cheapest way to avoid tuning a limit that was never responsible. + +Expect a binding limit to be **bidirectional**: truncating a stage wins on instances where its +effort was wasted and loses where it was load-bearing. A near-zero aggregate delta can hide +large per-instance swings in both directions, so judge by the split, never by the mean alone. + +## Log one parseable line per stage + +Emit the features the budget was derived from, the budget that came out, and what was actually +spent, on a single greppable line per stage. That makes a sweep regressable offline without +re-deriving anything from the solver, and it is what makes the attribution checks above +one-liners instead of re-runs. Include whether the limit bound and which of several limits +bound, since "the budget was 300 and it spent 302" is the fact that distinguishes a limit that +shaped the stage from one that stopped it. From ee2bf78714d7af88784b4cf8646fe3c11ecfb79d Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 23:38:23 +0200 Subject: [PATCH 07/20] Loosen presolve limits that cost reduction without saving needed time The cost policy capped Papilo rounds at 30 and clamped the probing badge above a 2e5 cost proxy. Measuring both directly shows neither limit was buying time that was actually scarce: - mzzv11 keeps reducing past 30 rounds (1962 rows against 1576) for 6s more, while triptim1 is bit-identical at 30 and unlimited. A round cap costs reduction where it binds and buys nothing where it does not, so rounds are now uncapped and the wall ceiling bounds the cost. - triptim1 (proxy 3.5e5) removes 669 rows at badge 1024 against 545 at 32, for 17s instead of 10s against a 60s ceiling. The threshold moves to 5e5 so the clamp only applies where the badge genuinely cannot fit, as on ns1760995 which needs 100s+ to make it pay. The structural probing ceiling is demoted to a backstop. Its cost proxy predicts throughput only to within ~700x, so no scale both bounds the wall time and leaves useful coverage: the value that kept every run under 60s also truncated the large instances to 1-9% of their candidates, which cost netdiversion and roi5alpha10n8 their solutions. Stopping at a measured wall target dominates it, so max_time_on_probing becomes the real bound and the ceiling only binds where the proxy is extreme. Verified at 600s: mzzv11 3.48 -> 0.26 error, triptim1 1.63 -> 0.00, netdiversion 41.70 -> 0.00, all better than baseline. The probing change is not yet settled on the noisier instances; configs 1, 2 and 3 span it for the full run. --- .../mip/heuristics_hyper_params.hpp | 14 +++-- cpp/src/math_optimization/solver_settings.cu | 2 +- .../diversity/diversity_config.hpp | 6 +- .../presolve/presolve_budget_policy.hpp | 56 +++++++++++-------- 4 files changed, 47 insertions(+), 31 deletions(-) diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index 5c0d6ebbc8..4c55b9e208 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -37,12 +37,14 @@ struct mip_heuristics_hyper_params_t { // differs between instances. f_t probe_host_overhead_work = 0.02; // charged per probed variable f_t probe_iter_work = 0.01; // charged per multi-probe propagation iteration - // Numerator of the work ceiling that keeps probing inside a wall-clock budget. Work units are - // reproducible but their cost is not: realised throughput spans 1.9 to 689 units/s across the - // 240-instance benchmark, so a coverage target alone lets slow instances run for the whole solve. - // Dividing this by the cost proxy (nnz + n_cand * avg_col_len) bounds the wall time instead; at - // 1.5e8 the slowest of 660 measured probing runs finishes in 44s, none above 60s. - f_t probing_work_time_scale = 1.5e8; + // Numerator of the probing work ceiling, divided by the cost proxy (nnz + n_cand * avg_col_len). + // This is a reproducible backstop for pathological instances, not the wall-clock bound: the proxy + // predicts realised throughput only to within ~700x, so a scale tight enough to bound time + // starves everything else. At 1.5e8 no run exceeded 44s but the large instances got 1-9% of their + // candidates probed, which cost several of them their solution; buying even 14% coverage there + // needs 1e9, which puts 22 runs past 60s. max_time_on_probing bounds the time instead, measured + // rather than predicted, and this stays loose enough to bind only where the proxy is extreme. + f_t probing_work_time_scale = 3.0e9; f_t root_lp_time_ratio = 0.1; // fraction of total time for root LP f_t root_lp_max_time = 15.0; // hard cap on root LP seconds diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index 370cd0b3f6..d23d59fdc9 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -109,7 +109,7 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings {CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT, &mip_settings.heuristic_params.cuopt_presolve_work_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(30.0), "probing-cache budget in work units (manual policy, or ceiling for derived policies)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK, &mip_settings.heuristic_params.probe_host_overhead_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.02), "work units charged per probed variable (host overhead)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK, &mip_settings.heuristic_params.probe_iter_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.01), "work units charged per multi-probe propagation iteration"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(1.5e8), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0e9), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, diff --git a/cpp/src/mip_heuristics/diversity/diversity_config.hpp b/cpp/src/mip_heuristics/diversity/diversity_config.hpp index ec6998c464..d5ec9e215b 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_config.hpp +++ b/cpp/src/mip_heuristics/diversity/diversity_config.hpp @@ -13,7 +13,11 @@ namespace cuopt::mathematical_optimization::mip { struct diversity_config_t { double time_ratio_of_probing_cache = 0.1; - double max_time_on_probing = 60.0; + // The bound that actually holds on probing. Measuring beats predicting here: stopping at a wall + // target leaves ~80% of instances probing every candidate, where a structural work ceiling tight + // enough to guarantee the same time truncates the large ones to a few percent. 45s keeps margin + // under the 60s presolve budget. + double max_time_on_probing = 45.0; int max_var_diff = 256; double default_time_limit = 10.; int initial_island_size = 3; diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index 724a6c37ef..e524da411e 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -69,16 +69,17 @@ struct presolve_config_t { }; inline constexpr presolve_config_t presolve_configs[] = { - // Old Papilo rule, so the probing ceiling is measured on its own. - {presolve_budget_policy_t::fixed, 1.00, 1.5e8}, // 0 - // New Papilo cost rule against the ceiling and the two coverage fractions. - {presolve_budget_policy_t::cost, 1.00, 1.5e8}, // 1 - {presolve_budget_policy_t::cost, 0.25, 1.5e8}, // 2 - {presolve_budget_policy_t::cost, 0.05, 1.5e8}, // 3 - // Half the ceiling: 29s worst case rather than 44s, to see whether the margin costs quality. - {presolve_budget_policy_t::cost, 1.00, 7.5e7}, // 4 - // Coverage only, no ceiling -- the control that says whether the ceiling earns its keep. - {presolve_budget_policy_t::cost, 0.05, 0.0}, // 5 + // Old Papilo rule, so the Papilo change can be read against the probing change. + {presolve_budget_policy_t::fixed, 1.00, 3.0e9}, // 0 + // The defaults. + {presolve_budget_policy_t::cost, 1.00, 3.0e9}, // 1 + // The tight ceiling that over-truncated: bounded time, but 1-9% coverage on large instances. + {presolve_budget_policy_t::cost, 1.00, 1.5e8}, // 2 + // No ceiling at all -- isolates how much the backstop still contributes over the wall target. + {presolve_budget_policy_t::cost, 1.00, 0.0}, // 3 + // Coverage fractions, to check the wall target is not leaving cheap probing on the table. + {presolve_budget_policy_t::cost, 0.25, 3.0e9}, // 4 + {presolve_budget_policy_t::cost, 0.05, 3.0e9}, // 5 }; inline constexpr int n_presolve_configs = 6; @@ -249,21 +250,30 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(feat.avg_col_len(), 1.0); - b.papilo_max_rounds = 30; - b.papilo_max_badgesize = papilo_probe_cost > 2.0e5 ? 32 : -1; + b.papilo_max_rounds = -1; + b.papilo_max_badgesize = papilo_probe_cost > 5.0e5 ? 32 : -1; probe_fraction = 1.0; b.probing_step_size = 128; break; From 3a2ba290671df69c5819deb0140f72d00a1b9c9c Mon Sep 17 00:00:00 2001 From: akif Date: Thu, 30 Jul 2026 23:39:04 +0200 Subject: [PATCH 08/20] Correct the bind-detection advice in the stage budgets reference The reference recommended reading hit_tlim to decide whether a round cap binds. That field reports the time limit only, so it stays clear on a stage the round cap stopped early -- following the advice exonerated a cap that was in fact costing reduction. Replace it with an A/B against the cap removed, and note that a binding round cap is better deleted than tuned since the wall ceiling bounds cost directly. Also add a noise-band check as the first attribution step. Using the spread between repeats of the same build as a per-instance threshold removed 18 of 26 apparent regressions on a 240-instance run. --- .../references/stage_budgets.md | 25 ++++++++++++++++--- 1 file changed, 21 insertions(+), 4 deletions(-) diff --git a/skills/cuopt-developer/references/stage_budgets.md b/skills/cuopt-developer/references/stage_budgets.md index 7858d81496..be0444998b 100644 --- a/skills/cuopt-developer/references/stage_budgets.md +++ b/skills/cuopt-developer/references/stage_budgets.md @@ -34,10 +34,20 @@ get less than they could use. ## Verify which knob costs, rather than the one that looks like it should Iteration or round counts read like cost limits and often are not: a stage that converges -before its round cap is not bounded by it at all. Check whether the cap ever binds -(`rounds_used < rounds_cap`, `hit_tlim=0`) before tuning it. The knob that drives cost is -usually the one controlling work *per* round — batch or candidate-set size — and it may sit -behind a third-party default. +before its round cap is not bounded by it at all. The knob that drives cost is usually the one +controlling work *per* round — batch or candidate-set size — and it may sit behind a +third-party default. + +Establish whether a cap binds by **running the stage with it removed and diffing the output**, +not by reading a flag. Flags answer narrower questions than they appear to: a `hit_tlim` / +`timed_out` field reports the *time* limit only, so it stays clear on a stage the round cap +stopped early, and reading it as "no limit bound" silently exonerates the wrong knob. Unless a +counter records the cap's own trigger, the A/B is the only sound check. + +When a round cap does bind, prefer deleting it over tuning it. It is a poor cost limit in both +directions — where it binds it truncates a stage mid-convergence and costs output, and where it +does not it saves nothing — while the wall ceiling bounds cost directly and only fires when +time is genuinely short. A knob can also be free on some instances and a genuine quality lever on others. Sweep it and record both cost and the stage's own output (reduced dimensions, cuts kept) — identical output @@ -55,6 +65,13 @@ gap): that pairing means the stage starved the root relaxation. Run these checks in order; each can end the investigation: +0. **Is the regression larger than the instance's own noise?** Use the spread between repeats + of the *same* build as a per-instance noise band and drop everything inside it. This is not + a formality: on one 240-instance run it removed 18 of 26 apparent regressions, and four of + the remaining eight turned out to have a bit-identical reduced problem, leaving four real + ones. Tuning against noise is worse than not tuning, because the change looks justified. + Instances with wide bands cannot be settled by single runs at all — carry them to the full + set as competing configs rather than picking a winner from one measurement. 1. **Did the limit bind?** Split instances by the `budget_exhausted` flag and compare the mean error delta per group. If the regressed set is not the bound set, the limit is not the cause. From 4d0c6a53e6731c97348c19b483b8a6618128fc0d Mon Sep 17 00:00:00 2001 From: akif Date: Fri, 31 Jul 2026 11:22:37 +0200 Subject: [PATCH 09/20] Bound the exempt Papilo badge and cut probing coverage to a quarter Two findings from the six-config sweep on 3a2ba2906. Exempting an instance from the badge clamp handed Papilo an unbounded badge (ncols/2), which on wide problems is enormous: rail01 ran at badge 58763 and Papilo 12.1s -> 49.2s, tbfp-network at 36373 and 4.8s -> 47.7s, neos-5114902-kasavu at 355082 and 7.3s -> 44.6s. Eight instances lost 12-43s of their presolve budget that way. Exempt now means 1024, which still keeps the reduction that badge 1024 bought on triptim1 and mzzv11. This also corrects the previous attribution. Configs 0 and 1 moved the round cap and the badge together, so the 231-instance net of +0.02 read as "the Papilo change does nothing" when the badge was doing the damage and the rounds were doing the good -- mzzv11 is 0.07 error uncapped against 3.41 capped, and none of the eight blowups came from rounds. Rounds stay uncapped and each sweep point now moves one knob. Probing coverage drops to a quarter, which beat full coverage by 0.615 mean error over 231 instances. Probing takes its time from branch and bound, and the wins (satellites2-40 100 -> 31.6, brazil3, rail01, physiciansched3-3) outweigh the losses on 30n20b8 and netdiversion. Not addressed here: square47 runs 2250s against a 600s limit in every config while presolve costs it under 40s, exploring 4-10 nodes where the baseline managed 2 in 600s, and its error is identical at 6.25 throughout. Cheap presolve just lets it start more giant node LPs with no deadline check inside one. That belongs in the node-LP time check. --- .../presolve/presolve_budget_policy.hpp | 76 ++++++++++++------- 1 file changed, 49 insertions(+), 27 deletions(-) diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index e524da411e..8d34f6c2a2 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -52,34 +52,38 @@ inline const char* presolve_budget_policy_name(int policy) } // Benchmark points selected by CUOPT_CONFIG_ID, overriding the policy hyper-parameter so one build -// covers the whole sweep. The previous 240-instance run showed that neither stage was actually -// bounded: presolve ran 243-403s (Papilo, ns1760995) and up to 598s (probing, six instances) out of -// a 600s solve, and the instances that lost their dual bound are the ones that spent it there. Both -// stages now carry a cost model instead of a coverage target, so these points span the two models -// against the previously measured coverage fractions. +// covers the whole sweep. // -// `work_time_scale` of 0 disables the probing ceiling, which separates "probe a fixed fraction" -// from "probe until a wall-clock proxy is exhausted"; the fractions bracket the crossover the -// earlier run pointed at, where truncation won on bab6 (0.5% of candidates), square41 (4.3%) and -// square47 (2.9%) and lost on 30n20b8 (21.6%) and physiciansched3-3 (9.1%). +// Per-instance error on this set is far too noisy to tune against: across two baseline repeats +// rail01 spans 20.6 to 35.8 and satellites2-40 spans 31.6 to 68.4, so any single instance can move +// tens of points for free. Only the mean over the whole set, and Papilo/probing wall time (which is +// stable), carry signal -- so each point below moves exactly one knob against config 0. struct presolve_config_t { presolve_budget_policy_t papilo_rule; + int papilo_rounds; // -1 uncapped + int badge_clamp; // badge at or above the cost threshold; <=0 applies no clamp + int badge_exempt; // badge below it; -1 uncapped, i.e. Papilo's ncols/2 double probe_fraction; double work_time_scale; }; +// One knob moves per point, all against config 0, so each arm is a paired A/B on the same build. +// The previous sweep could not attribute its own result: config 0 changed the round cap and the +// badge together, and the 231-instance net came out at +0.02, hiding that the badge was doing all +// the damage while the rounds were doing all the good. inline constexpr presolve_config_t presolve_configs[] = { - // Old Papilo rule, so the Papilo change can be read against the probing change. - {presolve_budget_policy_t::fixed, 1.00, 3.0e9}, // 0 - // The defaults. - {presolve_budget_policy_t::cost, 1.00, 3.0e9}, // 1 - // The tight ceiling that over-truncated: bounded time, but 1-9% coverage on large instances. - {presolve_budget_policy_t::cost, 1.00, 1.5e8}, // 2 - // No ceiling at all -- isolates how much the backstop still contributes over the wall target. - {presolve_budget_policy_t::cost, 1.00, 0.0}, // 3 - // Coverage fractions, to check the wall target is not leaving cheap probing on the table. - {presolve_budget_policy_t::cost, 0.25, 3.0e9}, // 4 - {presolve_budget_policy_t::cost, 0.05, 3.0e9}, // 5 + // The defaults, and the reference every other point is read against. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, // 0 + // Round cap restored. Uncapping is what took mzzv11 from 3.41 error to 0.07, and none of the + // Papilo blowups were rounds, but that was measured with the badge confounded -- so re-test it. + {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, // 1 + // Coverage either side of 0.25, which beat full coverage by 0.615 mean error over 231 instances. + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, // 2 + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.10, 3.0e9}, // 3 + // No clamp: is 32 on the expensive instances still worth it once the badge is bounded at 1024? + {presolve_budget_policy_t::cost, -1, 0, 1024, 0.25, 3.0e9}, // 4 + // No probing ceiling, leaving max_time_on_probing as the only bound. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 0.0}, // 5 }; inline constexpr int n_presolve_configs = 6; @@ -267,15 +271,26 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t 49.2s, tbfp-network to + // 36373 and 4.8s -> 47.7s, neos-5114902-kasavu to 355082 and 7.3s -> 44.6s. Eight instances + // regressed that way on one run, all of them from the badge rather than the round cap. 1024 is + // large enough to keep the reduction that mattered on triptim1 and mzzv11. + // + // Rounds are uncapped. A round cap looks like a cost limit but is not one: triptim1 is // bit-identical at 30 rounds and unlimited, while mzzv11 keeps reducing past 30 (1962 rows - // against 1576) for 6s more, so capping only cost reduction. The wall ceiling bounds the cost. + // against 1576, and 0.07 error against 3.41) for 6s more. None of the Papilo blowups above came + // from rounds, so there is nothing for the cap to save; the wall ceiling bounds the cost. case presolve_budget_policy_t::cost: { const double papilo_probe_cost = n_bin * std::max(feat.avg_col_len(), 1.0); b.papilo_max_rounds = -1; - b.papilo_max_badgesize = papilo_probe_cost > 5.0e5 ? 32 : -1; - probe_fraction = 1.0; - b.probing_step_size = 128; + b.papilo_max_badgesize = papilo_probe_cost > 5.0e5 ? 32 : 1024; + // Probing takes its time from branch and bound, and full coverage is not worth what it + // costs there: a quarter beat it by 0.615 mean error over 231 instances, winning big on + // satellites2-40 (100 -> 31.6), brazil3 and physiciansched3-3 against smaller losses. + probe_fraction = 0.25; + b.probing_step_size = 128; break; } @@ -313,8 +328,15 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t= 0) { - probe_fraction = presolve_configs[config].probe_fraction; - work_time_scale = presolve_configs[config].work_time_scale; + const auto& c = presolve_configs[config]; + probe_fraction = c.probe_fraction; + work_time_scale = c.work_time_scale; + // Re-derive the badge from the config's own pair so the sweep can separate "clamp the expensive + // instances" from "bound every badge", which the previous table conflated. + const double papilo_probe_cost = n_bin * std::max(feat.avg_col_len(), 1.0); + b.papilo_max_rounds = c.papilo_rounds; + b.papilo_max_badgesize = + (c.badge_clamp > 0 && papilo_probe_cost > 5.0e5) ? c.badge_clamp : c.badge_exempt; } probe_fraction *= work_scale; b.intended_probe_fraction = std::min(probe_fraction, 1.0); From 9ede4315e0be1f3a61b36bfb08ae95edb95a43e5 Mon Sep 17 00:00:00 2001 From: akif Date: Fri, 31 Jul 2026 12:28:31 +0200 Subject: [PATCH 10/20] Keep probing inside the presolve share of the solve The probing cap was min(max_time_on_probing, global_timer.remaining_time()), dropping the time_limit argument that main bounds it with. That argument is the presolve share of the solve (presolve_time_ratio, a tenth by default), so without it probing was bounded only by its own 45s ceiling and a short solve could spend all of itself there: at --time-limit 10 it could take the whole 10s where it should get 1s. At 600s the ceiling was tighter than main's 60s either way, which is why the benchmarks did not show it. Measured on supportcase6, probing now gets 0.63s at a 10s limit and 5.61s at 60s, and still reaches the 45s ceiling at 600s. --- .../mip_heuristics/diversity/diversity_manager.cu | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index a7e98e8773..b98c39df2a 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -322,10 +322,15 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ // shapes how much probing happens but cannot bound how long that takes: realised throughput // spans 12 to 450 work units per second across the benchmark, so at the slow end even a quarter // of the candidates ran 593 of the 600 available seconds and left the root LP with no dual - // bound. max_time_on_probing is therefore a hard ceiling rather than a target, and - // remaining_time only stops it reaching past the end of the solve. - f_t time_for_probing_cache = - std::min((f_t)diversity_config.max_time_on_probing, (f_t)global_timer.remaining_time()); + // bound. max_time_on_probing is therefore a hard ceiling rather than a target. + // + // time_limit is the presolve share of the solve (presolve_time_ratio, so a tenth by default) + // and has to stay in the minimum: dropping it left probing bounded only by its own 45s ceiling, + // so a short solve would spend all of itself here -- at --time-limit 10 probing could take the + // whole 10s where it should get 1s. remaining_time only stops it reaching past the end of the + // solve. + f_t time_for_probing_cache = std::min( + {(f_t)diversity_config.max_time_on_probing, time_limit, (f_t)global_timer.remaining_time()}); timer_t probing_timer{time_for_probing_cache}; const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem From 182b4a61b2f1ccd4bece4dc55de5d8ad942190c3 Mon Sep 17 00:00:00 2001 From: akif Date: Fri, 31 Jul 2026 12:57:33 +0200 Subject: [PATCH 11/20] Pair the sweep ids so each setting gets two repeats The harness runs every config id once and cannot repeat one, so the three settings worth measuring are duplicated across adjacent id pairs. In opportunistic mode a duplicated pair does not repeat itself; it measures the run-to-run spread, which is what the previous sweep lacked. The effects here are worth roughly 0.6 mean error while two baseline repeats differ by 0.53, so one run per setting cannot separate signal from noise. 0, 1 the defaults 2, 3 full coverage, against the quarter that beat it by 0.615 4, 5 round cap restored, now that the badge no longer confounds it Dropped from the previous table: the no-clamp arm, which square47 already settles with a bit-identical reduced problem at badge 32 against 1024, and the 0.10 coverage and no-ceiling arms, which only matter once the default holds up. --- .../presolve/presolve_budget_policy.hpp | 36 ++++++++++--------- 1 file changed, 20 insertions(+), 16 deletions(-) diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index 8d34f6c2a2..af3262b291 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -67,23 +67,27 @@ struct presolve_config_t { double work_time_scale; }; -// One knob moves per point, all against config 0, so each arm is a paired A/B on the same build. -// The previous sweep could not attribute its own result: config 0 changed the round cap and the -// badge together, and the 231-instance net came out at +0.02, hiding that the badge was doing all -// the damage while the rounds were doing all the good. +// Three settings, each duplicated across an adjacent pair of ids, so a harness that runs every id +// once still produces two repeats of each. Opportunistic mode is timing-dependent, so a duplicated +// pair does not repeat itself -- it measures the run-to-run spread, which is the thing the previous +// sweep lacked. The effects here are worth roughly 0.6 mean error while two baseline repeats differ +// by 0.53, so one run per setting cannot separate them. +// +// Only one knob moves per pair, against the pair at 0/1. The previous sweep could not attribute its +// own result: it changed the round cap and the badge together, and the 231-instance net came out at +// +0.02, hiding that the badge was doing all the damage while the rounds were doing all the good. inline constexpr presolve_config_t presolve_configs[] = { - // The defaults, and the reference every other point is read against. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, // 0 - // Round cap restored. Uncapping is what took mzzv11 from 3.41 error to 0.07, and none of the - // Papilo blowups were rounds, but that was measured with the badge confounded -- so re-test it. - {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, // 1 - // Coverage either side of 0.25, which beat full coverage by 0.615 mean error over 231 instances. - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, // 2 - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.10, 3.0e9}, // 3 - // No clamp: is 32 on the expensive instances still worth it once the badge is bounded at 1024? - {presolve_budget_policy_t::cost, -1, 0, 1024, 0.25, 3.0e9}, // 4 - // No probing ceiling, leaving max_time_on_probing as the only bound. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 0.0}, // 5 + // 0, 1: the defaults, and the reference the other pairs are read against. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, + // 2, 3: full coverage. A quarter beat it by 0.615 mean error over 231 instances, on one run -- + // this is the largest effect measured here, so it gets a second look before it is trusted. + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, + // 4, 5: round cap restored. Uncapping rests on mzzv11 alone (0.07 error against 3.41) and was + // measured with the badge confounded, which makes it the weakest-evidenced of the defaults. + {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, + {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, }; inline constexpr int n_presolve_configs = 6; From e3261f10fa31ee9c692ab11b17e51a35ab1103e8 Mon Sep 17 00:00:00 2001 From: akif Date: Fri, 31 Jul 2026 17:47:22 +0200 Subject: [PATCH 12/20] Sweep the probing work ceiling with the wall cap removed The 45s cap on probing was a stand-in for a bound the work budget was meant to provide, so let a config remove it and size the ceiling instead. The presolve share of the solve and the global timer stay: they bound probing against the rest of the solve rather than tuning it. Three arms set the coverage target to 1, which carries no target at all, so the work limit is exactly the structural ceiling and the scale-to-time curve is readable without the fraction masking which of the two bound an instance. Two more pair the ceiling with the usual quarter coverage as default candidates. Papilo is pinned across all of them so the probing bound is the only thing moving. --- .../diversity/diversity_manager.cu | 11 +++- .../presolve/presolve_budget_policy.hpp | 63 ++++++++++++------- 2 files changed, 51 insertions(+), 23 deletions(-) diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index b98c39df2a..285af235aa 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -329,8 +329,15 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ // so a short solve would spend all of itself here -- at --time-limit 10 probing could take the // whole 10s where it should get 1s. remaining_time only stops it reaching past the end of the // solve. - f_t time_for_probing_cache = std::min( - {(f_t)diversity_config.max_time_on_probing, time_limit, (f_t)global_timer.remaining_time()}); + // + // A sweep config replaces the dedicated ceiling, so an arm can test whether the work budget + // bounds probing on its own. time_limit and remaining_time are not part of that experiment: + // they bound probing against the rest of the solve rather than tuning it, so they stay. + const f_t probing_wall_cap = probing_budget.config_id >= 0 + ? (f_t)probing_budget.probing_wall_limit + : (f_t)diversity_config.max_time_on_probing; + f_t time_for_probing_cache = + std::min({probing_wall_cap, time_limit, (f_t)global_timer.remaining_time()}); timer_t probing_timer{time_for_probing_cache}; const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index af3262b291..9a9c57a13f 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -65,29 +65,43 @@ struct presolve_config_t { int badge_exempt; // badge below it; -1 uncapped, i.e. Papilo's ncols/2 double probe_fraction; double work_time_scale; + // Dedicated wall cap on probing in seconds; <=0 removes it and leaves the work budget as the only + // bound the policy imposes. + double probing_wall_limit; }; -// Three settings, each duplicated across an adjacent pair of ids, so a harness that runs every id -// once still produces two repeats of each. Opportunistic mode is timing-dependent, so a duplicated -// pair does not repeat itself -- it measures the run-to-run spread, which is the thing the previous -// sweep lacked. The effects here are worth roughly 0.6 mean error while two baseline repeats differ -// by 0.53, so one run per setting cannot separate them. +// This sweep asks one question: can the work budget bound probing on its own, with the dedicated +// wall cap removed? Everything except the probing bound is pinned at the current default, so any +// difference is attributable to it. // -// Only one knob moves per pair, against the pair at 0/1. The previous sweep could not attribute its -// own result: it changed the round cap and the badge together, and the 231-instance net came out at -// +0.02, hiding that the badge was doing all the damage while the rounds were doing all the good. +// Removing the cap does not leave probing unbounded. The presolve share of the solve stays in the +// minimum -- a tenth of the limit, so 60s at the benchmark's 600s -- and that is the number these +// arms are measured against: an arm whose ceiling is too loose shows up as probing sitting at 60s +// with a truncated candidate set, not as the 553s it would take unbounded. The count of instances +// that reach it is the result, and 0 means the ceiling did the job. +// +// The scale is the knob under test because worst-case probing time is close to linear in it: 1.5e8 +// held the worst run to 44.7s, 1e9 reached 297s and 3e9 reached 553s. Tight is not free, though -- +// at 1.5e8 the large instances got 1-9% of their candidates probed and several lost their solution +// -- so the arms bracket the range rather than assuming the safe end is usable. +// +// Unlike the previous table these are six distinct points rather than three duplicated pairs. What +// is being measured here is probing wall time, whether the bound was reached, and coverage, all of +// which are stable across repeats; it was per-instance *error* that needed the pairing, and error +// is only a sanity check on the mean here. inline constexpr presolve_config_t presolve_configs[] = { - // 0, 1: the defaults, and the reference the other pairs are read against. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9}, - // 2, 3: full coverage. A quarter beat it by 0.615 mean error over 231 instances, on one run -- - // this is the largest effect measured here, so it gets a second look before it is trusted. - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 3.0e9}, - // 4, 5: round cap restored. Uncapping rests on mzzv11 alone (0.07 error against 3.41) and was - // measured with the badge confounded, which makes it the weakest-evidenced of the defaults. - {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, - {presolve_budget_policy_t::cost, 30, 32, 1024, 0.25, 3.0e9}, + // 0: control. The current default, wall cap included, and the reference for the others. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9, 45.0}, + // 1, 2, 3: the ceiling as the sole work bound. A fraction at or above 1 carries no coverage + // target, so probing_work_limit is exactly the ceiling and the scale-to-time curve is visible + // without the fraction masking which of the two bound a given instance. + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 0.0}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 4.0e8, 0.0}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.0e9, 0.0}, + // 4, 5: the candidates for a default. Coverage bounds the cheap instances and the ceiling bounds + // the expensive ones, which is the split the wall cap was standing in for. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 4.0e8, 0.0}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 1.0e9, 0.0}, }; inline constexpr int n_presolve_configs = 6; @@ -152,6 +166,10 @@ struct presolve_budget_t { // what set probing_work_limit. Logged so a run can be attributed to one or the other offline. double probing_work_ceiling{std::numeric_limits::infinity()}; bool probing_ceiling_binding{false}; + // Dedicated wall cap on probing, infinite when the work budget is meant to be the only bound. + // The presolve share of the solve and the global timer still apply on top of this; they are + // correctness bounds rather than tuning knobs and no config removes them. + double probing_wall_limit{std::numeric_limits::infinity()}; // The policy that actually ran and the config that selected it, so a log line is attributable // even when CUOPT_CONFIG_ID overrode the hyper-parameter. int policy{1}; @@ -341,6 +359,8 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t 0 && papilo_probe_cost > 5.0e5) ? c.badge_clamp : c.badge_exempt; + b.probing_wall_limit = + c.probing_wall_limit > 0.0 ? c.probing_wall_limit : std::numeric_limits::infinity(); } probe_fraction *= work_scale; b.intended_probe_fraction = std::min(probe_fraction, 1.0); @@ -369,7 +389,7 @@ inline void log_presolve_budget(const char* stage, "PRESOLVE_BUDGET stage=%s config=%d policy=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f " "nbin=%.0f arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " "rounds=%d badge=%d work=%.3f step=%d intended_probe_frac=%.4f work_ceiling=%.3f " - "ceiling_binding=%d", + "ceiling_binding=%d wall=%.1f", stage, b.config_id, presolve_budget_policy_name(b.policy), @@ -390,7 +410,8 @@ inline void log_presolve_budget(const char* stage, b.probing_step_size, b.intended_probe_fraction, b.probing_work_ceiling, - (int)b.probing_ceiling_binding); + (int)b.probing_ceiling_binding, + b.probing_wall_limit); } } // namespace cuopt::mathematical_optimization::mip From 4f71e4aed36debcb67c74b933ca4a6196816639e Mon Sep 17 00:00:00 2001 From: akif Date: Mon, 3 Aug 2026 15:01:38 +0200 Subject: [PATCH 13/20] Make the probing work ceiling the bound and demote the wall cap Adopts the sweep's tightest arm as the default. Over 240 instances a scale of 1.5e8 stopped every probing run before the wall cap could fire, worst case 42.8s, where 4e8 and 1e9 still needed the wall on 2 and 8 instances. It also spends the least time doing it -- 582s of probing across the set against 1136s at 4e8 -- and solves that overran their own limit fell from 17 to 7. max_time_on_probing goes to 120s as a backstop for an instance the cost proxy badly mispredicts. Leaving it tight would keep it shaping the common case, which is what the work ceiling now does. Drops the coverage fraction so the ceiling is the sole work bound. A fraction cuts every instance by the same proportion regardless of cost, which both obscured attribution and penalised cheap instances; probing now runs to completion on those while the ceiling truncates by cost. Mean error is unchanged within noise (11.78 against 12.01, with 0.53 between repeats of one build), so this was decided on bounding behaviour and time. --- .../mip/heuristics_hyper_params.hpp | 25 +++++--- cpp/src/math_optimization/solver_settings.cu | 2 +- .../diversity/diversity_config.hpp | 13 +++-- .../diversity/diversity_manager.cu | 29 +++++----- .../presolve/presolve_budget_policy.hpp | 57 +++++++++---------- 5 files changed, 68 insertions(+), 58 deletions(-) diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index 4c55b9e208..91f6830a4e 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -38,13 +38,24 @@ struct mip_heuristics_hyper_params_t { f_t probe_host_overhead_work = 0.02; // charged per probed variable f_t probe_iter_work = 0.01; // charged per multi-probe propagation iteration // Numerator of the probing work ceiling, divided by the cost proxy (nnz + n_cand * avg_col_len). - // This is a reproducible backstop for pathological instances, not the wall-clock bound: the proxy - // predicts realised throughput only to within ~700x, so a scale tight enough to bound time - // starves everything else. At 1.5e8 no run exceeded 44s but the large instances got 1-9% of their - // candidates probed, which cost several of them their solution; buying even 14% coverage there - // needs 1e9, which puts 22 runs past 60s. max_time_on_probing bounds the time instead, measured - // rather than predicted, and this stays loose enough to bind only where the proxy is extreme. - f_t probing_work_time_scale = 3.0e9; + // This is the bound on probing, not a backstop: measured over 240 instances it stopped every run + // before the wall cap could fire, with a worst case of 42.8s, while the arms at 4e8 and 1e9 + // needed the wall on 2 and 8 instances respectively. It also spends the least time to do it -- + // 582s of probing across the set against 1136s at 4e8 -- and instances whose solve overran the + // limit fell from 17 to 7. + // + // Tight enough to bound time is also tight enough to truncate: 34 instances end below 5% coverage + // here against 20 under the older loose scale. That trade is deliberate and measured neutral on + // solution quality (11.78 mean error against 12.01, inside the 0.53 run-to-run noise), because + // probing takes its time from branch and bound. + // + // A larger scale is not the way to buy coverage back. The proxy predicts throughput only to + // within ~340x and one instance (nw04) pins the scale under every reshaping tried, including + // refitting the exponent to the measured nnz^0.65. Giving wide-row problems (avg row length above + // ~100) their own scale is worth 2.3x on this one, and beyond that the residual is not explained + // by any structural feature -- it needs throughput measured during probing rather than predicted + // from the problem. + f_t probing_work_time_scale = 1.5e8; f_t root_lp_time_ratio = 0.1; // fraction of total time for root LP f_t root_lp_max_time = 15.0; // hard cap on root LP seconds diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index d23d59fdc9..370cd0b3f6 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -109,7 +109,7 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings {CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT, &mip_settings.heuristic_params.cuopt_presolve_work_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(30.0), "probing-cache budget in work units (manual policy, or ceiling for derived policies)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK, &mip_settings.heuristic_params.probe_host_overhead_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.02), "work units charged per probed variable (host overhead)"}, {CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK, &mip_settings.heuristic_params.probe_iter_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.01), "work units charged per multi-probe propagation iteration"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0e9), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(1.5e8), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, diff --git a/cpp/src/mip_heuristics/diversity/diversity_config.hpp b/cpp/src/mip_heuristics/diversity/diversity_config.hpp index d5ec9e215b..92d14bba40 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_config.hpp +++ b/cpp/src/mip_heuristics/diversity/diversity_config.hpp @@ -13,11 +13,14 @@ namespace cuopt::mathematical_optimization::mip { struct diversity_config_t { double time_ratio_of_probing_cache = 0.1; - // The bound that actually holds on probing. Measuring beats predicting here: stopping at a wall - // target leaves ~80% of instances probing every candidate, where a structural work ceiling tight - // enough to guarantee the same time truncates the large ones to a few percent. 45s keeps margin - // under the 60s presolve budget. - double max_time_on_probing = 45.0; + // Backstop only. probing_work_time_scale is what bounds probing -- across 240 instances it + // stopped every run before this could fire, worst case 42.8s -- so this exists for an instance + // whose cost the proxy misses badly, not for routine use. Kept loose deliberately: tightening it + // would start shaping the common case again, which is what the work ceiling replaced. + // + // At the default presolve share (a tenth of the solve) this cannot bind below a 1200s limit, + // since the share itself is the smaller bound before then. + double max_time_on_probing = 120.0; int max_var_diff = 256; double default_time_limit = 10.; int initial_island_size = 3; diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 285af235aa..162451b633 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -300,7 +300,7 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ const auto& hp = context.settings.heuristic_params; const bool deterministic = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC; const auto probing_features = probing_presolve_features(*problem_ptr); - const auto probing_budget = evaluate_presolve_budget(hp, probing_features); + auto probing_budget = evaluate_presolve_budget(hp, probing_features); bool run_probing_cache = !fj_only_run; // Under the legacy policy the probing cache carries no work budget, so in deterministic mode // nothing would bound it: the wall clock is infinite there. Keep it off, which is also what makes @@ -317,27 +317,26 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ run_probing_cache = false; } if (run_probing_cache) { - log_presolve_budget("PROBING", probing_features, probing_budget); // Run probing cache before trivial presolve to discover variable implications. The work budget - // shapes how much probing happens but cannot bound how long that takes: realised throughput - // spans 12 to 450 work units per second across the benchmark, so at the slow end even a quarter - // of the candidates ran 593 of the 600 available seconds and left the root LP with no dual - // bound. max_time_on_probing is therefore a hard ceiling rather than a target. + // is what bounds this: at probing_work_time_scale it stopped every one of 240 instances before + // max_time_on_probing could fire, so the wall is a backstop for a badly mispredicted instance + // rather than the working limit it used to be. // // time_limit is the presolve share of the solve (presolve_time_ratio, so a tenth by default) - // and has to stay in the minimum: dropping it left probing bounded only by its own 45s ceiling, + // and has to stay in the minimum: without it probing would be bounded only by its own wall cap, // so a short solve would spend all of itself here -- at --time-limit 10 probing could take the // whole 10s where it should get 1s. remaining_time only stops it reaching past the end of the // solve. // - // A sweep config replaces the dedicated ceiling, so an arm can test whether the work budget - // bounds probing on its own. time_limit and remaining_time are not part of that experiment: - // they bound probing against the rest of the solve rather than tuning it, so they stay. - const f_t probing_wall_cap = probing_budget.config_id >= 0 - ? (f_t)probing_budget.probing_wall_limit - : (f_t)diversity_config.max_time_on_probing; - f_t time_for_probing_cache = - std::min({probing_wall_cap, time_limit, (f_t)global_timer.remaining_time()}); + // A sweep config replaces the wall cap so an arm can be measured with the work budget as the + // only bound. time_limit and remaining_time are not part of that experiment: they bound probing + // against the rest of the solve rather than tuning it, so they stay. + if (probing_budget.config_id < 0) { + probing_budget.probing_wall_limit = diversity_config.max_time_on_probing; + } + log_presolve_budget("PROBING", probing_features, probing_budget); + f_t time_for_probing_cache = std::min( + {(f_t)probing_budget.probing_wall_limit, time_limit, (f_t)global_timer.remaining_time()}); timer_t probing_timer{time_for_probing_cache}; const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index 9a9c57a13f..ae6fc33ef9 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -70,38 +70,32 @@ struct presolve_config_t { double probing_wall_limit; }; -// This sweep asks one question: can the work budget bound probing on its own, with the dedicated -// wall cap removed? Everything except the probing bound is pinned at the current default, so any -// difference is attributable to it. +// The scale ladder that settled the probing bound, kept so the result can be re-measured. Slot 0 is +// the adopted default; the rest are the arms it was chosen over. Only the probing bound moves. // -// Removing the cap does not leave probing unbounded. The presolve share of the solve stays in the -// minimum -- a tenth of the limit, so 60s at the benchmark's 600s -- and that is the number these -// arms are measured against: an arm whose ceiling is too loose shows up as probing sitting at 60s -// with a truncated candidate set, not as the 553s it would take unbounded. The count of instances -// that reach it is the result, and 0 means the ceiling did the job. +// Result over 240 instances: the number of instances where probing hit a wall bound was 0 at 1.5e8, +// 2 at 4e8 and 8 at 1e9, and total probing time across the set was 582s, 1136s and 1867s. Solves +// that overran their own limit followed the same ordering (7, 27, 24) against 17 for the +// wall-capped arm. Mean error spanned 11.48 to 12.51 across all six arms while two repeats of one +// build differ by 0.53, so error could not rank them and did not decide this. // -// The scale is the knob under test because worst-case probing time is close to linear in it: 1.5e8 -// held the worst run to 44.7s, 1e9 reached 297s and 3e9 reached 553s. Tight is not free, though -- -// at 1.5e8 the large instances got 1-9% of their candidates probed and several lost their solution -// -- so the arms bracket the range rather than assuming the safe end is usable. -// -// Unlike the previous table these are six distinct points rather than three duplicated pairs. What -// is being measured here is probing wall time, whether the bound was reached, and coverage, all of -// which are stable across repeats; it was per-instance *error* that needed the pairing, and error -// is only a sanity check on the mean here. +// Note what these arms do *not* test. The presolve share of the solve stays in the minimum, so an +// arm whose ceiling is too loose shows up as probing sitting at that share with a truncated +// candidate set rather than running unbounded. Only 1.5e8 is bounded by its own ceiling; 4e8 and +// 1e9 are still leaning on a wall, which is the real reason they were rejected. inline constexpr presolve_config_t presolve_configs[] = { - // 0: control. The current default, wall cap included, and the reference for the others. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9, 45.0}, - // 1, 2, 3: the ceiling as the sole work bound. A fraction at or above 1 carries no coverage - // target, so probing_work_limit is exactly the ceiling and the scale-to-time curve is visible - // without the fraction masking which of the two bound a given instance. - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 0.0}, + // 0, 1: the adopted default. No coverage target, so the ceiling is the sole work bound, and the + // wall left loose enough that it is a backstop rather than the thing shaping probing. + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 120.0}, + {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 120.0}, + // 2, 3: looser scales. Rejected for needing the wall on 2 and 8 instances, and for spending 2-3x + // the probing time to do it. {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 4.0e8, 0.0}, {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.0e9, 0.0}, - // 4, 5: the candidates for a default. Coverage bounds the cheap instances and the ceiling bounds - // the expensive ones, which is the split the wall cap was standing in for. + // 4: a quarter coverage alongside the ceiling, the shape used before this sweep. {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 4.0e8, 0.0}, - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 1.0e9, 0.0}, + // 5: the previous default, wall-capped and with a scale loose enough never to bind. + {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9, 45.0}, }; inline constexpr int n_presolve_configs = 6; @@ -308,10 +302,13 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t(feat.avg_col_len(), 1.0); b.papilo_max_rounds = -1; b.papilo_max_badgesize = papilo_probe_cost > 5.0e5 ? 32 : 1024; - // Probing takes its time from branch and bound, and full coverage is not worth what it - // costs there: a quarter beat it by 0.615 mean error over 231 instances, winning big on - // satellites2-40 (100 -> 31.6), brazil3 and physiciansched3-3 against smaller losses. - probe_fraction = 0.25; + // No coverage target: the work ceiling is the only bound, which is what let the wall cap go. + // A coverage fraction cuts every instance by the same proportion whether or not it is the + // expensive one, so pairing it with the ceiling only obscured which of the two stopped a + // given run. Dropping it also stops penalising the cheap instances -- probing runs over a + // second on 102 instances here against 81 before, while runs over ten seconds fell from 34 to + // 11, since the ceiling truncates by cost rather than uniformly. + probe_fraction = 1.0; b.probing_step_size = 128; break; } From 0b77cd173d2c9124f5dc65cf0f9f836a27f49d97 Mon Sep 17 00:00:00 2001 From: akif Date: Mon, 3 Aug 2026 17:49:59 +0200 Subject: [PATCH 14/20] Report a too-dense ADAT as a capacity failure instead of crashing rd-rplusc-21 was the 239-of-240 instance, dying in the concurrent barrier LP with "Requested size overflows device_uvector storage" and taking the whole MIP solve with it, since the error escaped run_mip and suppressed the result row. It was never a large allocation, it was a negative one. The instance is 125899 x 622 with four columns each touching ~119000 rows, so A*A^T is 90.5% dense at 1.44e10 nonzeros. Dense-column elimination does its job -- the barrier finds 8, and dropping 8 brings A*A^T to 1.0e9 -- but the residual system still needs 2311527529 nonzeros, clearing INT32_MAX by 7.6%. cuSPARSE reports that count as int64_t while resize_to_nnz takes i_t, so it narrowed to -1983439767 and reached RMM as ~1.8e19. That 7.6% margin is why the failure looked nondeterministic. Small differences in what presolve leaves behind put the count either side of the limit, and below it the same instance instead asks for ~28 GB and fails as an already-handled out_of_memory. Guarding the narrowing reports the real count and cause, and the barrier's handler now catches std::bad_alloc rather than only rmm::out_of_memory so any allocation failure there degrades to NUMERICAL_ISSUES, which is how the barrier already treats every other shape it cannot hold. The instance is still hard -- 0 nodes explored in 600s -- but it now terminates normally and reports a result rather than vanishing from the set. Not addressed here: the augmented-system fallback cannot fire for this shape, needing estimated_nz_AAT > 1e10 and max_row_nz > 5000 while max row length is 100 and m^2 caps the estimate at 2.9e9. --- cpp/src/barrier/barrier.cu | 5 ++++- cpp/src/barrier/sparse_matrix_kernels.cuh | 17 +++++++++++++++++ 2 files changed, 21 insertions(+), 1 deletion(-) diff --git a/cpp/src/barrier/barrier.cu b/cpp/src/barrier/barrier.cu index 784f6c0901..502267ed4c 100644 --- a/cpp/src/barrier/barrier.cu +++ b/cpp/src/barrier/barrier.cu @@ -4439,7 +4439,10 @@ lp_status_t barrier_solver_t::solve(f_t start_time, lp_solution_t #include +#include + +#include +#include + namespace cuopt::mathematical_optimization::barrier { template @@ -140,6 +145,18 @@ void multiply_kernels(raft::handle_t const* handle, int64_t ADAT_num_rows, ADAT_num_cols, ADAT_nnz1; RAFT_CUSPARSE_TRY( cusparseSpMatGetSize(cusparse_data.matADAT_descr, &ADAT_num_rows, &ADAT_num_cols, &ADAT_nnz1)); + // cuSPARSE sizes the product in 64 bits, but the CSR arrays are indexed by i_t. A tall problem + // with a few near-full columns leaves ADAT dense enough to pass that range even after dense + // columns are eliminated, and narrowing would hand RMM a negative count that resurfaces as an + // unrelated "size overflows device_uvector storage" from deep inside the allocator. Report the + // real cause instead, as the capacity failure the caller already knows how to fall back from. + if (ADAT_nnz1 > std::numeric_limits::max()) { + throw rmm::out_of_memory( + "ADAT needs " + std::to_string(ADAT_nnz1) + " nonzeros over " + + std::to_string(ADAT_num_rows) + " rows, past the " + + std::to_string(std::numeric_limits::max()) + + " its index type can address: the normal equations are too dense for this problem"); + } ADAT.resize_to_nnz(ADAT_nnz1, handle->get_stream()); thrust::fill(rmm::exec_policy(handle->get_stream()), ADAT.x.begin(), ADAT.x.end(), 0.0); From 837a37e94ade8ba6bc04318bbc34a2d2b7c2b1fa Mon Sep 17 00:00:00 2001 From: akif Date: Tue, 4 Aug 2026 11:09:07 +0200 Subject: [PATCH 15/20] Keep the measured presolve rule and drop the sweep scaffolding The policy menu existed to choose between competing budget rules on a benchmark; that question is settled, so the enum, the six-arm config table and the CUOPT_CONFIG_ID selector go, leaving the one rule that was adopted. Default behaviour is unchanged by construction: rounds uncapped, badge 32 or 1024 on the same n_bin * avg_col_len > 5e5 threshold, work limit 1.5e8 / (nnz + n_cand * avg_col_len), step 128, wall 120s -- the wall now set unconditionally rather than gated on the config id. The work model moves from tunable hyper-parameters to constants beside the scale they calibrate, since the ceiling is only meaningful against the costs the probing loop actually charges. Removes the settings for the policy, the step size, both work weights, the work-time scale, and cuopt_presolve_work_limit along with the coverage multiplier it drove, so the ceiling is the only work bound. presolve_max_rounds and papilo_probing_max_badgesize stay as overrides, now defaulting to -1 meaning "derive from the problem", where an explicit 0 removes the cap. The per-stage Papilo telemetry drops to DEBUG now that it is no longer being regressed offline, and papilo_presolve_features moves next to the Papilo code it feeds rather than sitting in solve.cu. --- .../utils/presolve_budget_sweep.py | 310 -------------- .../mathematical_optimization/constants.h | 10 - .../mip/heuristics_hyper_params.hpp | 37 +- cpp/src/math_optimization/solver_settings.cu | 10 +- .../diversity/diversity_manager.cu | 8 +- .../presolve/presolve_budget_policy.hpp | 387 ++++-------------- .../mip_heuristics/presolve/probing_cache.cu | 9 +- .../presolve/third_party_presolve.cpp | 30 +- .../presolve/third_party_presolve.hpp | 6 + cpp/src/mip_heuristics/solve.cu | 32 +- 10 files changed, 117 insertions(+), 722 deletions(-) delete mode 100644 benchmarks/linear_programming/utils/presolve_budget_sweep.py diff --git a/benchmarks/linear_programming/utils/presolve_budget_sweep.py b/benchmarks/linear_programming/utils/presolve_budget_sweep.py deleted file mode 100644 index 0c5e82298a..0000000000 --- a/benchmarks/linear_programming/utils/presolve_budget_sweep.py +++ /dev/null @@ -1,310 +0,0 @@ -# SPDX-FileCopyrightText: Copyright (c) 2026, NVIDIA CORPORATION & AFFILIATES. All rights reserved. -# SPDX-License-Identifier: Apache-2.0 -r"""Sweep the MIP presolve budget policies over a set of instances. - -Runs cuopt_cli once per (instance, policy), scrapes the PRESOLVE_* log lines plus the final result -line, and writes one CSV row per run. The CSV carries both the structural features each budget was -derived from and what the budget actually spent, which is what a fit of -"budget <- problem dimensions and structure" needs. - -Sweeping whole policies:: - - python presolve_budget_sweep.py --cli ./cpp/build/cuopt_cli \ - --dataset-dir datasets/mip/miplib2017 --time-limit 300 --policies 0 2 \ - --out /tmp/presolve_sweep.csv --log-dir /tmp/presolve_sweep_logs - -Sweeping one knob instead gives the mapping curve from a wall-clock limit onto a -round / badge / work-unit budget:: - - python presolve_budget_sweep.py --grid-param mip_hyper_heuristic_presolve_max_rounds \ - --grid-values 1 2 3 5 10 20 30 50 --time-limit 2 --out /tmp/rounds_map.csv -""" - -import argparse -import csv -import os -import re -import subprocess -import sys -import time - -# Chosen to span the structural axes the policies key on: nnz, average row length, and binary -# fraction. bab2 / supportcase6 / square41 are the instances whose presolve ran away unbounded and -# motivated the budgets in the first place, so they are the regression cases. -DEFAULT_INSTANCES = [ - "bab2", - "supportcase6", - "square41", - "air05", - "30n20b8", - "gen-ip054", - "nw04", - "rail507", - "seymour", - "mzzv11", - "roll3000", - "ns1208400", - "glass4", - "timtab1", - "enlight_hard", - "sp97ar", -] - -POLICY_NAMES = { - 0: "legacy", - 1: "fixed", - 2: "size", - 3: "density", - 4: "binary", - 5: "combined", - 6: "manual", -} - -KV_RE = re.compile(r"(\w+)=(-?[\w.+-]+)") -EXPLORED_RE = re.compile( - r"Explored (\d+) nodes \((\d+) simplex iterations\) in ([\d.]+)s" -) -OBJ_RE = re.compile( - r"Best objective ([-\d.eE+]+), best bound ([-\d.eE+]+), gap ([-\d.eE+]+|inf)%" -) - -STATUS_MARKERS = [ - ("Optimal solution found", "Optimal"), - ("Time limit reached", "TimeLimit"), - ("Work limit reached", "WorkLimit"), - ("Problem is infeasible", "Infeasible"), - ("Problem is unbounded", "Unbounded"), - ("No solution found", "NoSolution"), - # A claim of integer infeasibility on an instance known to be feasible means a reduction was - # unsound, so it must be distinguishable from simply not finding a solution in time. - ("Problem has no integer feasible solution", "NoIntegerFeasible"), -] - - -def parse_kv(line, prefix): - """Pull every key=value token that follows `prefix` on `line`.""" - idx = line.find(prefix) - if idx < 0: - return {} - return dict(KV_RE.findall(line[idx + len(prefix) :])) - - -def parse_log(text): - row = {} - for line in text.splitlines(): - if "PRESOLVE_BUDGET stage=PAPILO" in line: - for k, v in parse_kv(line, "PRESOLVE_BUDGET").items(): - row[f"papilo_{k}"] = v - elif "PRESOLVE_BUDGET stage=PROBING" in line: - for k, v in parse_kv(line, "PRESOLVE_BUDGET").items(): - row[f"probing_{k}"] = v - elif "PRESOLVE_PAPILO_REDUCED" in line: - for k, v in parse_kv(line, "PRESOLVE_PAPILO_REDUCED").items(): - row[f"reduced_{k}"] = v - elif "PRESOLVE_PAPILO wall=" in line: - for k, v in parse_kv(line, "PRESOLVE_PAPILO").items(): - row["papilo_wall" if k == "wall" else f"papilo_{k}"] = v - elif "PRESOLVE_PROBING_WALL" in line: - row["probing_wall"] = parse_kv(line, "PRESOLVE_PROBING_WALL").get( - "wall" - ) - elif "PRESOLVE_PROBING probes=" in line: - for k, v in parse_kv(line, "PRESOLVE_PROBING").items(): - row[f"spent_{k}"] = v - elif "Probing-cache step disabled" in line: - row["probing_disabled"] = "1" - - m = EXPLORED_RE.search(line) - if m: - row["nodes"], row["simplex_iters"], row["solve_wall"] = m.groups() - m = OBJ_RE.search(line) - if m: - row["objective"], row["bound"], row["gap_pct"] = m.groups() - for marker, status in STATUS_MARKERS: - if marker in line: - row["status"] = status - return row - - -def as_text(stream): - """Decode a stream, since subprocess hands back bytes on the timeout path even with text=True.""" - if stream is None: - return "" - if isinstance(stream, bytes): - return stream.decode("utf-8", "replace") - return stream - - -def run_one(args, instance, policy, config_path, grid_value=None): - with open(config_path, "w") as fh: - fh.write(f"mip_hyper_heuristic_presolve_budget_policy = {policy}\n") - if grid_value is not None: - fh.write(f"{args.grid_param} = {grid_value}\n") - for extra in args.param: - fh.write(extra.replace(":", " = ", 1) + "\n") - - cmd = [ - args.cli, - os.path.join(args.dataset_dir, instance + ".mps"), - "--time-limit", - str(args.time_limit), - "--params-file", - config_path, - ] - if args.determinism: - cmd += ["--mip-determinism-mode", "1"] - - t0 = time.time() - timed_out = False - try: - proc = subprocess.run( - cmd, - capture_output=True, - text=True, - timeout=args.timeout, - check=False, - ) - output = proc.stdout + proc.stderr - returncode = proc.returncode - except subprocess.TimeoutExpired as exc: - output = as_text(exc.stdout) + as_text(exc.stderr) - returncode = -1 - timed_out = True - wall = time.time() - t0 - - row = parse_log(output) - row.update( - instance=instance, - policy=policy, - policy_name=POLICY_NAMES.get(policy, str(policy)), - harness_wall=f"{wall:.2f}", - harness_timeout=int(timed_out), - returncode=returncode, - ) - if grid_value is not None: - row["grid_param"] = args.grid_param - row["grid_value"] = grid_value - if timed_out: - row.setdefault("status", "HarnessTimeout") - - if args.log_dir: - os.makedirs(args.log_dir, exist_ok=True) - with open( - os.path.join(args.log_dir, f"{instance}.p{policy}.log"), "w" - ) as fh: - fh.write(output) - return row - - -def main(): - ap = argparse.ArgumentParser() - ap.add_argument("--cli", default="./cpp/build/cuopt_cli") - ap.add_argument("--dataset-dir", default="datasets/mip/miplib2017") - ap.add_argument("--instances", nargs="*", default=DEFAULT_INSTANCES) - ap.add_argument( - "--policies", nargs="*", type=int, default=[0, 1, 2, 3, 4, 5] - ) - ap.add_argument("--time-limit", type=float, default=30.0) - ap.add_argument( - "--timeout", - type=float, - default=300.0, - help="hard wall cap per run; the legacy policy leaves presolve unbounded", - ) - ap.add_argument("--determinism", action="store_true") - ap.add_argument( - "--param", - action="append", - default=[], - help="extra config entry as key:value, repeatable", - ) - ap.add_argument( - "--grid-param", - default="", - help="sweep this single hyper-parameter under the manual policy instead of sweeping " - "policies; this is what yields the wall-limit -> rounds / work-unit mapping curve", - ) - ap.add_argument("--grid-values", nargs="*", default=[]) - ap.add_argument("--out", default="presolve_sweep.csv") - ap.add_argument("--log-dir", default="") - args = ap.parse_args() - - config_path = "/tmp/presolve_budget_sweep.config" - rows = [] - # A grid sweeps one knob under the manual policy; otherwise the variant axis is the policy. - if args.grid_param: - variants = [(6, v) for v in args.grid_values] - else: - variants = [(p, None) for p in args.policies] - - total = len(args.instances) * len(variants) - done = 0 - for instance in args.instances: - path = os.path.join(args.dataset_dir, instance + ".mps") - if not os.path.exists(path): - print(f"SKIP missing {path}", flush=True) - continue - for policy, grid_value in variants: - done += 1 - # A sweep is long enough that losing all of it to one bad run is the worst outcome. - try: - row = run_one(args, instance, policy, config_path, grid_value) - except Exception as exc: # noqa: BLE001 - row = { - "instance": instance, - "policy": policy, - "policy_name": POLICY_NAMES.get(policy, str(policy)), - "grid_value": grid_value, - "status": "HarnessError", - "harness_wall": "0", - "harness_error": repr(exc), - } - rows.append(row) - variant = ( - f"{args.grid_param.split('_')[-1]}={grid_value}" - if grid_value is not None - else f"p{policy} {row.get('policy_name', '')}" - ) - print( - f"[{done}/{total}] {instance:16s} {variant:14s}" - f" status={row.get('status', '?'):18s}" - f" papilo_wall={row.get('papilo_wall', '-'):>8s}" - f" red_vars={row.get('reduced_nvars', '-'):>8s}" - f" red_nnz={row.get('reduced_nnz', '-'):>9s}" - f" probing_wall={row.get('probing_wall', '-'):>8s}" - f" probes={row.get('spent_probes', '-'):>7s}" - f" work={row.get('spent_work', '-'):>9s}", - flush=True, - ) - # Written incrementally so a long sweep is inspectable while it runs. - write_csv(args.out, rows) - - write_csv(args.out, rows) - print(f"\nwrote {len(rows)} rows to {args.out}") - - -def write_csv(path, rows): - if not rows: - return - fields = [] - for row in rows: - for key in row: - if key not in fields: - fields.append(key) - lead = [ - "instance", - "policy", - "policy_name", - "grid_param", - "grid_value", - "status", - ] - fields = lead + [f for f in fields if f not in lead] - with open(path, "w", newline="") as fh: - writer = csv.DictWriter(fh, fieldnames=fields, extrasaction="ignore") - writer.writeheader() - writer.writerows(rows) - - -if __name__ == "__main__": - sys.exit(main()) diff --git a/cpp/include/cuopt/mathematical_optimization/constants.h b/cpp/include/cuopt/mathematical_optimization/constants.h index 22917c3431..00ef468946 100644 --- a/cpp/include/cuopt/mathematical_optimization/constants.h +++ b/cpp/include/cuopt/mathematical_optimization/constants.h @@ -100,19 +100,9 @@ #define CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS "mip_hyper_heuristic_num_cpufj_threads" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO "mip_hyper_heuristic_presolve_time_ratio" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME "mip_hyper_heuristic_presolve_max_time" -#define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY \ - "mip_hyper_heuristic_presolve_budget_policy" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS "mip_hyper_heuristic_presolve_max_rounds" #define CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE \ "mip_hyper_heuristic_papilo_probing_max_badgesize" -#define CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT \ - "mip_hyper_heuristic_cuopt_presolve_work_limit" -#define CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE "mip_hyper_heuristic_probing_step_size" -#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK \ - "mip_hyper_heuristic_probe_host_overhead_work" -#define CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK "mip_hyper_heuristic_probe_iter_work" -#define CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE \ - "mip_hyper_heuristic_probing_work_time_scale" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO "mip_hyper_heuristic_root_lp_time_ratio" #define CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME "mip_hyper_heuristic_root_lp_max_time" #define CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT "mip_hyper_heuristic_rins_time_limit" diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index 91f6830a4e..94ed27f492 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -24,38 +24,11 @@ struct mip_heuristics_hyper_params_t { f_t presolve_time_ratio = 0.1; // fraction of total time for presolve f_t presolve_max_time = 60.0; // hard cap on presolve seconds - // Presolve budgeting. presolve_budget_policy selects how the four knobs below are derived from - // the problem's dimensions and structure (see presolve_budget_policy.hpp); the values here are - // the defaults the policy starts from and the literal values used by the `manual` policy. - i_t presolve_budget_policy = 7; // presolve_budget_policy_t (cost) - i_t presolve_max_rounds = 30; // Papilo presolve rounds cap (<=0 = Papilo default) - i_t papilo_probing_max_badgesize = 1024; // ceiling on Papilo's probing.minbadgesize - f_t cuopt_presolve_work_limit = 30.0; // probing-cache budget, work units - i_t probing_step_size = 512; // probed vars between work-budget checks - // Weights of the probing-cache work model. Work units measure probing effort reproducibly; they - // are not an estimate of elapsed time, so the effort-per-second they correspond to legitimately - // differs between instances. - f_t probe_host_overhead_work = 0.02; // charged per probed variable - f_t probe_iter_work = 0.01; // charged per multi-probe propagation iteration - // Numerator of the probing work ceiling, divided by the cost proxy (nnz + n_cand * avg_col_len). - // This is the bound on probing, not a backstop: measured over 240 instances it stopped every run - // before the wall cap could fire, with a worst case of 42.8s, while the arms at 4e8 and 1e9 - // needed the wall on 2 and 8 instances respectively. It also spends the least time to do it -- - // 582s of probing across the set against 1136s at 4e8 -- and instances whose solve overran the - // limit fell from 17 to 7. - // - // Tight enough to bound time is also tight enough to truncate: 34 instances end below 5% coverage - // here against 20 under the older loose scale. That trade is deliberate and measured neutral on - // solution quality (11.78 mean error against 12.01, inside the 0.53 run-to-run noise), because - // probing takes its time from branch and bound. - // - // A larger scale is not the way to buy coverage back. The proxy predicts throughput only to - // within ~340x and one instance (nw04) pins the scale under every reshaping tried, including - // refitting the exponent to the measured nnz^0.65. Giving wide-row problems (avg row length above - // ~100) their own scale is worth 2.3x on this one, and beyond that the residual is not explained - // by any structural feature -- it needs throughput measured during probing rather than predicted - // from the problem. - f_t probing_work_time_scale = 1.5e8; + // Presolve budgeting. Both are derived from the problem's dimensions and structure by default + // (see presolve_budget_policy.hpp); a negative value asks for that rule, and any other value + // overrides it, where <=0 removes the cap entirely. + i_t presolve_max_rounds = -1; // Papilo presolve rounds cap + i_t papilo_probing_max_badgesize = -1; // ceiling on Papilo's probing.minbadgesize f_t root_lp_time_ratio = 0.1; // fraction of total time for root LP f_t root_lp_max_time = 15.0; // hard cap on root LP seconds diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index ef55735d65..181949b2d7 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -106,10 +106,6 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO, &mip_settings.heuristic_params.presolve_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for presolve"}, {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME, &mip_settings.heuristic_params.presolve_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(60.0), "hard cap on presolve seconds"}, - {CUOPT_MIP_HYPER_HEURISTIC_CUOPT_PRESOLVE_WORK_LIMIT, &mip_settings.heuristic_params.cuopt_presolve_work_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(30.0), "probing-cache budget in work units (manual policy, or ceiling for derived policies)"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBE_HOST_OVERHEAD_WORK, &mip_settings.heuristic_params.probe_host_overhead_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.02), "work units charged per probed variable (host overhead)"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBE_ITER_WORK, &mip_settings.heuristic_params.probe_iter_work, f_t(0.0), std::numeric_limits::infinity(), f_t(0.01), "work units charged per multi-probe propagation iteration"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBING_WORK_TIME_SCALE, &mip_settings.heuristic_params.probing_work_time_scale, f_t(0.0), std::numeric_limits::infinity(), f_t(1.5e8), "numerator of the probing work ceiling; divided by nnz + n_cand * avg_col_len (0 disables)"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, @@ -172,10 +168,8 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") {CUOPT_MIP_HYPER_HEURISTIC_POPULATION_SIZE, &mip_settings.heuristic_params.population_size, 1, std::numeric_limits::max(), 32, "max solutions in pool"}, {CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS, &mip_settings.heuristic_params.num_cpufj_threads, 0, std::numeric_limits::max(), 8, "parallel CPU FJ climbers"}, - {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_BUDGET_POLICY, &mip_settings.heuristic_params.presolve_budget_policy, 0, 7, 7, "how presolve budgets are derived: 0=legacy 1=fixed 2=size 3=density 4=binary 5=combined 6=manual 7=cost"}, - {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS, &mip_settings.heuristic_params.presolve_max_rounds, -1, std::numeric_limits::max(), 30, "Papilo presolve rounds cap (<=0 keeps Papilo default)"}, - {CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE, &mip_settings.heuristic_params.papilo_probing_max_badgesize, -1, std::numeric_limits::max(), 1024, "ceiling on Papilo probing.minbadgesize (<=0 leaves it uncapped)"}, - {CUOPT_MIP_HYPER_HEURISTIC_PROBING_STEP_SIZE, &mip_settings.heuristic_params.probing_step_size, 1, std::numeric_limits::max(), 512, "probed variables between probing-cache work-budget checks"}, + {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS, &mip_settings.heuristic_params.presolve_max_rounds, -1, std::numeric_limits::max(), -1, "Papilo presolve rounds cap (<0 derives it from the problem, 0 keeps Papilo default)"}, + {CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE, &mip_settings.heuristic_params.papilo_probing_max_badgesize, -1, std::numeric_limits::max(), -1, "ceiling on Papilo probing.minbadgesize (<0 derives it from the problem, 0 leaves it uncapped)"}, {CUOPT_MIP_HYPER_HEURISTIC_STAGNATION_TRIGGER, &mip_settings.heuristic_params.stagnation_trigger, 1, std::numeric_limits::max(), 3, "FP loops w/o improvement before recombination"}, {CUOPT_MIP_HYPER_HEURISTIC_MAX_ITERS_WITHOUT_IMPROVEMENT, &mip_settings.heuristic_params.max_iterations_without_improvement, 1, std::numeric_limits::max(), 8, "diversity step depth after stagnation"}, {CUOPT_MIP_HYPER_HEURISTIC_N_OF_MINIMUMS_FOR_EXIT, &mip_settings.heuristic_params.n_of_minimums_for_exit, 1, std::numeric_limits::max(), 7000, "FJ baseline local-minima exit threshold"}, diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 162451b633..64699a3102 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -327,13 +327,7 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ // so a short solve would spend all of itself here -- at --time-limit 10 probing could take the // whole 10s where it should get 1s. remaining_time only stops it reaching past the end of the // solve. - // - // A sweep config replaces the wall cap so an arm can be measured with the work budget as the - // only bound. time_limit and remaining_time are not part of that experiment: they bound probing - // against the rest of the solve rather than tuning it, so they stay. - if (probing_budget.config_id < 0) { - probing_budget.probing_wall_limit = diversity_config.max_time_on_probing; - } + probing_budget.probing_wall_limit = diversity_config.max_time_on_probing; log_presolve_budget("PROBING", probing_features, probing_budget); f_t time_for_probing_cache = std::min( {(f_t)probing_budget.probing_wall_limit, time_limit, (f_t)global_timer.remaining_time()}); diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index ae6fc33ef9..e8fea9565f 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -12,116 +12,51 @@ #include #include -#include -#include -#include #include -#include namespace cuopt::mathematical_optimization::mip { -// Competing hypotheses for how much presolve effort a problem deserves. Each policy maps the -// structural features below onto the same four budgets, so a sweep over policies on a fixed -// instance set isolates the effect of the mapping itself. `legacy` reproduces pre-budget behaviour -// and is the baseline to measure against; `manual` reads the hyper-parameters verbatim so a -// specific point can be pinned from the command line. -enum class presolve_budget_policy_t : int { - legacy = 0, - fixed = 1, - size = 2, - density = 3, - binary = 4, - combined = 5, - manual = 6, - cost = 7, -}; - -inline const char* presolve_budget_policy_name(int policy) -{ - switch (static_cast(policy)) { - case presolve_budget_policy_t::legacy: return "legacy"; - case presolve_budget_policy_t::fixed: return "fixed"; - case presolve_budget_policy_t::size: return "size"; - case presolve_budget_policy_t::density: return "density"; - case presolve_budget_policy_t::binary: return "binary"; - case presolve_budget_policy_t::combined: return "combined"; - case presolve_budget_policy_t::manual: return "manual"; - case presolve_budget_policy_t::cost: return "cost"; - default: return "unknown"; - } -} +// Work the probing loop charges per unit of effort. These are exact counts rather than timings, +// which is what makes the budget reproducible, and probing_work_time_scale below is calibrated +// against them -- changing either invalidates it. +inline constexpr double probing_probe_work = 0.02; // per probed variable, host overhead +inline constexpr double probing_iter_work = 0.01; // per multi-probe propagation iteration -// Benchmark points selected by CUOPT_CONFIG_ID, overriding the policy hyper-parameter so one build -// covers the whole sweep. +// Numerator of the probing work ceiling, divided by the cost proxy. This is the bound on probing, +// not a backstop: measured over 240 instances it stopped every run before the wall cap could fire, +// with a worst case of 44.7s, while looser scales of 4e8 and 1e9 needed the wall on 2 and 8 +// instances and spent 2-3x the total probing time to do it. // -// Per-instance error on this set is far too noisy to tune against: across two baseline repeats -// rail01 spans 20.6 to 35.8 and satellites2-40 spans 31.6 to 68.4, so any single instance can move -// tens of points for free. Only the mean over the whole set, and Papilo/probing wall time (which is -// stable), carry signal -- so each point below moves exactly one knob against config 0. -struct presolve_config_t { - presolve_budget_policy_t papilo_rule; - int papilo_rounds; // -1 uncapped - int badge_clamp; // badge at or above the cost threshold; <=0 applies no clamp - int badge_exempt; // badge below it; -1 uncapped, i.e. Papilo's ncols/2 - double probe_fraction; - double work_time_scale; - // Dedicated wall cap on probing in seconds; <=0 removes it and leaves the work budget as the only - // bound the policy imposes. - double probing_wall_limit; -}; - -// The scale ladder that settled the probing bound, kept so the result can be re-measured. Slot 0 is -// the adopted default; the rest are the arms it was chosen over. Only the probing bound moves. +// Tight enough to bound time is also tight enough to truncate, and that trade is deliberate: +// probing takes its time from branch and bound, and the truncation measured neutral on solution +// quality (11.78 mean error against 12.01, inside the 0.53 run-to-run noise). // -// Result over 240 instances: the number of instances where probing hit a wall bound was 0 at 1.5e8, -// 2 at 4e8 and 8 at 1e9, and total probing time across the set was 582s, 1136s and 1867s. Solves -// that overran their own limit followed the same ordering (7, 27, 24) against 17 for the -// wall-capped arm. Mean error spanned 11.48 to 12.51 across all six arms while two repeats of one -// build differ by 0.53, so error could not rank them and did not decide this. +// A larger scale is not the way to buy coverage back. The proxy predicts throughput only to within +// ~340x and one instance (nw04) pins the scale under every reshaping tried, including refitting the +// exponent to the measured nnz^0.65. Beyond that the residual is not explained by any structural +// feature; it needs throughput measured during probing rather than predicted from the problem. +inline constexpr double probing_work_time_scale = 1.5e8; + +// Probed variables between work-budget checks, i.e. the granularity at which the budget can be +// enforced. Work is only folded in at the step barrier, so too large a step runs unbudgeted. +inline constexpr int probing_budget_step_size = 128; + +// cuOpt forces probing.minbadgesize to ncols/2 to stop Papilo aborting probing on its own work +// budget, so the first badge overshoots that budget in a single pass and the badge cap is the only +// thing left bounding it. How much one badge costs scales with the nonzeros a probing sweep visits, +// n_bin * avg_col_len, so the clamp applies only where that cost is large enough to matter. // -// Note what these arms do *not* test. The presolve share of the solve stays in the minimum, so an -// arm whose ceiling is too loose shows up as probing sitting at that share with a truncated -// candidate set rather than running unbounded. Only 1.5e8 is bounded by its own ceiling; 4e8 and -// 1e9 are still leaning on a wall, which is the real reason they were rejected. -inline constexpr presolve_config_t presolve_configs[] = { - // 0, 1: the adopted default. No coverage target, so the ceiling is the sole work bound, and the - // wall left loose enough that it is a backstop rather than the thing shaping probing. - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 120.0}, - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.5e8, 120.0}, - // 2, 3: looser scales. Rejected for needing the wall on 2 and 8 instances, and for spending 2-3x - // the probing time to do it. - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 4.0e8, 0.0}, - {presolve_budget_policy_t::cost, -1, 32, 1024, 1.00, 1.0e9, 0.0}, - // 4: a quarter coverage alongside the ceiling, the shape used before this sweep. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 4.0e8, 0.0}, - // 5: the previous default, wall-capped and with a scale loose enough never to bind. - {presolve_budget_policy_t::cost, -1, 32, 1024, 0.25, 3.0e9, 45.0}, -}; -inline constexpr int n_presolve_configs = 6; - -// -1 when unset or out of range. Read once; the environment cannot change mid-run. -inline int presolve_config_id() -{ - static const int id = []() -> int { - const char* raw = std::getenv("CUOPT_CONFIG_ID"); - if (raw == nullptr) { return -1; } - try { - const int v = std::stoi(raw); - if (v < 0 || v >= n_presolve_configs) { - CUOPT_LOG_WARN("CUOPT_CONFIG_ID=%d is outside [0, %d); ignoring it for presolve budgets", - v, - n_presolve_configs); - return -1; - } - CUOPT_LOG_INFO("Using presolve budget config %d from CUOPT_CONFIG_ID", v); - return v; - } catch (const std::exception& e) { - CUOPT_LOG_WARN("Failed to parse CUOPT_CONFIG_ID: %s", e.what()); - return -1; - } - }(); - return id; -} +// Above this threshold the badge stops paying for itself: square47 (2.7e7) went from 40.6s at badge +// 1024 to 8.1s at 32 with a bit-identical reduced problem. Nearer the threshold the reduction is +// worth more than the time -- triptim1 (3.5e5) removes 669 rows at 1024 against 545 at 32, and +// clamping it regressed the instance -- and the same holds below it for mzzv11, 30n20b8 and air05. +// +// Being exempt still means a bounded badge. Leaving it uncapped hands Papilo ncols/2, which on wide +// problems is enormous and costs most of the presolve budget for nothing: rail01 went to badge +// 58763 and Papilo 12.1s -> 49.2s. 1024 keeps the reduction that mattered on triptim1 and mzzv11. +inline constexpr double papilo_badge_cost_threshold = 5.0e5; +inline constexpr int papilo_badge_clamped = 32; +inline constexpr int papilo_badge_exempt = 1024; // Dimensions plus cheap structural ratios. Both presolve stages populate this from whatever problem // representation they hold: Papilo from the original problem before any reduction, the cuOpt @@ -151,245 +86,68 @@ struct presolve_budget_t { // Probing-cache budget in work units: a reproducible count of probing effort, not a time // estimate. double probing_work_limit{std::numeric_limits::infinity()}; - // Probed variables per step, i.e. the granularity at which the budget can be enforced. - int probing_step_size{2048}; - // Coverage the policy asked for, kept only so the log can be compared against what was realised; - // the two diverge when an instance's probes are dearer than the average the budget assumed. - double intended_probe_fraction{1.0}; - // The work ceiling the cost proxy produced, and whether it rather than the coverage target is - // what set probing_work_limit. Logged so a run can be attributed to one or the other offline. - double probing_work_ceiling{std::numeric_limits::infinity()}; - bool probing_ceiling_binding{false}; - // Dedicated wall cap on probing, infinite when the work budget is meant to be the only bound. - // The presolve share of the solve and the global timer still apply on top of this; they are - // correctness bounds rather than tuning knobs and no config removes them. + int probing_step_size{probing_budget_step_size}; + // Dedicated wall cap on probing. The presolve share of the solve and the global timer still apply + // on top of this; they are correctness bounds rather than tuning knobs. double probing_wall_limit{std::numeric_limits::infinity()}; - // The policy that actually ran and the config that selected it, so a log line is attributable - // even when CUOPT_CONFIG_ID overrode the hyper-parameter. - int policy{1}; - int config_id{-1}; }; -namespace detail { - -inline double clamp_d(double v, double lo, double hi) { return std::min(std::max(v, lo), hi); } - -inline int clamp_i(double v, int lo, int hi) -{ - return static_cast(clamp_d(std::round(v), lo, hi)); -} - -} // namespace detail - +// Derives both presolve stages' budgets from the problem's dimensions and structure. Rounds and +// badge accept an explicit override from the hyper-parameters; a negative setting asks for the +// measured rule below, and any other value is passed through, where <=0 removes the cap. template presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t& hp, const presolve_features_t& feat) { - using detail::clamp_d; - using detail::clamp_i; - presolve_budget_t b{}; const double nnz = std::max(feat.nnz, 1.0); - const double arl = std::max(feat.avg_row_len(), 1.0); // Probing candidates are the integers of the problem the probing cache runs on. const double n_cand = std::max(feat.n_int, 1.0); const double n_bin = std::max(feat.n_bin, 1.0); - const double bf = feat.bin_frac(); + const double acl = std::max(feat.avg_col_len(), 1.0); + + // Rounds are uncapped. A round cap looks like a cost limit but is not one: triptim1 is + // bit-identical at 30 rounds and unlimited, while mzzv11 keeps reducing past 30 (1962 rows + // against 1576, and 0.07 error against 3.41) for 6s more. None of the Papilo blowups measured + // came from rounds, so there is nothing for the cap to save; the wall ceiling bounds the cost. + b.papilo_max_rounds = hp.presolve_max_rounds >= 0 ? (int)hp.presolve_max_rounds : -1; - // Probing cost is close to linear in the candidate count -- ~0.055 work units per candidate, - // measured across the 240-instance benchmark -- and the candidate set is the integers of the - // reduced problem. A constant budget therefore probes a small problem exhaustively and a large - // one barely at all: a flat 30 units covered 21% of 30n20b8 but 0.8% of netdiversion, which cost - // four instances their feasible solution. Policies state the fraction of candidates they want - // probed instead, converted below with the same cost model the probing loop charges. - constexpr double avg_iters_per_probe = 3.5; - const double per_candidate_work = - (double)hp.probe_host_overhead_work + avg_iters_per_probe * (double)hp.probe_iter_work; - // Lets a single knob scale every derived policy without recompiling: 1.0 at the default of 30. - const double work_scale = static_cast(hp.cuopt_presolve_work_limit) / 30.0; - double probe_fraction = 1.0; - double work_time_scale = (double)hp.probing_work_time_scale; + const double papilo_probe_cost = n_bin * acl; + b.papilo_max_badgesize = + hp.papilo_probing_max_badgesize >= 0 + ? (int)hp.papilo_probing_max_badgesize + : (papilo_probe_cost > papilo_badge_cost_threshold ? papilo_badge_clamped + : papilo_badge_exempt); // Cost of one probing sweep: every propagation touches the rows of the probed column, so the work // a second buys falls off with problem size. nnz + n_cand * avg_col_len tracked that better than - // nnz, arl or n_cand alone over 660 measured probing runs, and dividing the scale by it bounds - // the wall time that a pure coverage target cannot: throughput ranged 1.9 to 689 work units per - // second, so the same budget was worth 360x more time on one instance than another. - const double probing_cost_proxy = nnz + n_cand * std::max(feat.avg_col_len(), 1.0); - - // A config id, when set, replaces both the policy and its coverage. Resolved here rather than at - // each call site so the Papilo stage and the probing stage cannot disagree about which point of - // the sweep is running. - const int config = presolve_config_id(); - const auto policy = config >= 0 ? presolve_configs[config].papilo_rule - : (presolve_budget_policy_t)hp.presolve_budget_policy; - b.policy = (int)policy; - b.config_id = config; - - switch (policy) { - case presolve_budget_policy_t::legacy: - b.papilo_max_rounds = -1; - b.papilo_max_badgesize = -1; - b.probing_work_limit = std::numeric_limits::infinity(); - b.probing_step_size = 2048; - return b; - - // Papilo-only arm: the round and badge caps measured clean against the baseline, so probing is - // left unbounded here to isolate their effect from the probing budget's. - case presolve_budget_policy_t::fixed: - b.papilo_max_rounds = 30; - b.papilo_max_badgesize = 1024; - probe_fraction = 1.0; - b.probing_step_size = 128; - break; - - case presolve_budget_policy_t::manual: - b.papilo_max_rounds = hp.presolve_max_rounds; - b.papilo_max_badgesize = hp.papilo_probing_max_badgesize; - b.probing_work_limit = hp.cuopt_presolve_work_limit; - b.probing_step_size = hp.probing_step_size; - return b; - - // Measured rule. probing.minbadgesize, not the round count, is what drives Papilo's cost, and - // on wide problems a large badge buys almost nothing: capping it at 32 left the reduced problem - // bit-identical on square41/supportcase6/nw04/rail507 while cutting presolve 2-26x. On narrow - // problems the opposite holds (mzzv11, 30n20b8, air05 all reduce measurably worse at 32), so - // the cap is applied by width only. Rounds are kept non-binding on wide problems -- everything - // measured saturates well before 50 -- and capped on narrow ones, where mzzv11 keeps growing. - // The probing fraction is the hypothesis under test here, not a measured value: truncated - // probing was what won on bab6 (0.5% of candidates), square41 (4.3%) and square47 (2.9%), and - // what lost on 30n20b8 (21.6%) and physiciansched3-3 (9.1%), so the benchmark has one point per - // instance and no curve. A quarter sits between those two clusters. - case presolve_budget_policy_t::size: { - const bool wide = feat.n_vars > 2.0e4; - b.papilo_max_rounds = wide ? 50 : 20; - b.papilo_max_badgesize = wide ? 32 : -1; - probe_fraction = 0.25; - b.probing_step_size = 128; - break; - } - - // cuOpt forces probing.minbadgesize to ncols/2 to stop Papilo aborting probing on its own work - // budget (working_limit = 2*nnz), so the first badge overshoots that budget in a single pass - // and the badge cap is the only thing left bounding it. How much one badge costs scales with - // the nonzeros a probing sweep visits, n_bin * avg_col_len, but that predicts *cost*, not - // whether the badge is worth paying for, so the cap only applies where the cost is large enough - // to matter. - // - // Above ~5e5 the badge stops paying for itself: square47 (2.7e7) went from 40.6s at badge 1024 - // to 8.1s at 32 with a bit-identical reduced problem. ns1760995 (1.8e6) is the one instance - // where the badge does buy reduction -- it needs 640+ to remove 120k rows instead of 226 -- but - // that costs 100-274s against a 60s ceiling, and capping it matches the baseline gap (1.5% vs - // 1.4%) while letting it run cost every config its dual bound entirely. - // - // Nearer the threshold the reduction is worth more than the time: triptim1 (3.5e5) removes 669 - // rows at badge 1024 against 545 at 32, for 17s instead of 10s, and clamping it regressed the - // instance. Below it the same holds for mzzv11, 30n20b8 and air05. - // - // Being exempt from the clamp still means a bounded badge, not an unbounded one. Leaving it at - // -1 hands Papilo ncols/2, which on wide problems is enormous and costs most of the presolve - // budget for nothing: rail01 went to badge 58763 and Papilo 12.1s -> 49.2s, tbfp-network to - // 36373 and 4.8s -> 47.7s, neos-5114902-kasavu to 355082 and 7.3s -> 44.6s. Eight instances - // regressed that way on one run, all of them from the badge rather than the round cap. 1024 is - // large enough to keep the reduction that mattered on triptim1 and mzzv11. - // - // Rounds are uncapped. A round cap looks like a cost limit but is not one: triptim1 is - // bit-identical at 30 rounds and unlimited, while mzzv11 keeps reducing past 30 (1962 rows - // against 1576, and 0.07 error against 3.41) for 6s more. None of the Papilo blowups above came - // from rounds, so there is nothing for the cap to save; the wall ceiling bounds the cost. - case presolve_budget_policy_t::cost: { - const double papilo_probe_cost = n_bin * std::max(feat.avg_col_len(), 1.0); - b.papilo_max_rounds = -1; - b.papilo_max_badgesize = papilo_probe_cost > 5.0e5 ? 32 : 1024; - // No coverage target: the work ceiling is the only bound, which is what let the wall cap go. - // A coverage fraction cuts every instance by the same proportion whether or not it is the - // expensive one, so pairing it with the ceiling only obscured which of the two stopped a - // given run. Dropping it also stops penalising the cheap instances -- probing runs over a - // second on 102 instances here against 81 before, while runs over ten seconds fell from 34 to - // 11, since the ceiling truncates by cost rather than uniformly. - probe_fraction = 1.0; - b.probing_step_size = 128; - break; - } - - // A single propagation sweep costs roughly one pass over each touched row, so long rows make - // every probe more expensive: buy fewer of them, and size the badge so one badge's working - // limit (~2*nnz in Papilo) stays roughly constant. - case presolve_budget_policy_t::density: - b.papilo_max_rounds = arl <= 10 ? 40 : arl <= 50 ? 20 : 8; - b.papilo_max_badgesize = clamp_i(2.0e6 / arl, 32, 4096); - probe_fraction = clamp_d(10.0 / arl, 0.05, 1.0); - b.probing_step_size = arl <= 50 ? 256 : 128; - break; - - // Probing and clique merging only pay off on binaries, so scale with how many there are and how - // much of the problem they make up. - case presolve_budget_policy_t::binary: - b.papilo_max_rounds = bf >= 0.9 ? 50 : bf >= 0.5 ? 30 : 15; - b.papilo_max_badgesize = clamp_i(std::max(n_bin / 2.0, 32.0), 32, 1024); - probe_fraction = clamp_d(bf, 0.05, 1.0); - b.probing_step_size = 128; - break; - - // Multiplicative over the three effects above, so no single feature can dominate the budget. - case presolve_budget_policy_t::combined: { - const double size_f = clamp_d(1.0e5 / nnz, 0.25, 2.0); - const double density_f = clamp_d(10.0 / arl, 0.25, 2.0); - const double binary_f = clamp_d(0.5 + bf, 0.5, 1.5); - const double factor = size_f * density_f * binary_f; - b.papilo_max_rounds = clamp_i(30.0 * factor, 5, 60); - b.papilo_max_badgesize = clamp_i(2.0e6 / arl, 32, 2048); - probe_fraction = clamp_d(0.25 * factor, 0.05, 1.0); - b.probing_step_size = clamp_i(128.0 * density_f, 64, 512); - break; - } - } - - if (config >= 0) { - const auto& c = presolve_configs[config]; - probe_fraction = c.probe_fraction; - work_time_scale = c.work_time_scale; - // Re-derive the badge from the config's own pair so the sweep can separate "clamp the expensive - // instances" from "bound every badge", which the previous table conflated. - const double papilo_probe_cost = n_bin * std::max(feat.avg_col_len(), 1.0); - b.papilo_max_rounds = c.papilo_rounds; - b.papilo_max_badgesize = - (c.badge_clamp > 0 && papilo_probe_cost > 5.0e5) ? c.badge_clamp : c.badge_exempt; - b.probing_wall_limit = - c.probing_wall_limit > 0.0 ? c.probing_wall_limit : std::numeric_limits::infinity(); - } - probe_fraction *= work_scale; - b.intended_probe_fraction = std::min(probe_fraction, 1.0); - - // A fraction at or above 1 means "probe everything", carried as no budget rather than as a large - // number: per_candidate_work is an average, so on an instance whose probes are dearer than - // average a finite budget sized for full coverage would still cut probing short. - const double coverage_work = probe_fraction >= 1.0 ? std::numeric_limits::infinity() - : probe_fraction * n_cand * per_candidate_work; - b.probing_work_ceiling = work_time_scale > 0.0 ? work_time_scale / probing_cost_proxy - : std::numeric_limits::infinity(); - // Whichever binds first. The ceiling is what keeps probing away from the wall cap it used to run - // into, and the fraction is what stops it probing more of a cheap problem than is worth probing. - b.probing_ceiling_binding = b.probing_work_ceiling < coverage_work; - b.probing_work_limit = std::min(coverage_work, b.probing_work_ceiling); + // nnz, avg_row_len or n_cand alone over 660 measured probing runs, and dividing the scale by it + // bounds the wall time that a pure coverage target cannot: throughput ranged 1.9 to 689 work + // units per second, so the same budget was worth 360x more time on one instance than another. + // + // The ceiling is the only work bound; there is no coverage target alongside it. A coverage + // fraction cuts every instance by the same proportion whether or not it is the expensive one, + // which both obscured which of the two stopped a run and penalised the cheap instances. Dropping + // it left probing running over a second on 102 instances against 81 before, while runs over ten + // seconds fell from 34 to 11, since the ceiling truncates by cost rather than uniformly. + const double probing_cost_proxy = nnz + n_cand * acl; + b.probing_work_limit = probing_work_time_scale / probing_cost_proxy; + b.probing_step_size = probing_budget_step_size; return b; } // One line per presolve stage carrying the features that went in and the budgets that came out, so -// a sweep can be regressed offline without re-deriving anything from the solver. +// a run can be regressed offline without re-deriving anything from the solver. inline void log_presolve_budget(const char* stage, const presolve_features_t& f, const presolve_budget_t& b) { CUOPT_LOG_INFO( - "PRESOLVE_BUDGET stage=%s config=%d policy=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f " + "PRESOLVE_BUDGET stage=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f " "nbin=%.0f arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " - "rounds=%d badge=%d work=%.3f step=%d intended_probe_frac=%.4f work_ceiling=%.3f " - "ceiling_binding=%d wall=%.1f", + "rounds=%d badge=%d work=%.3f step=%d wall=%.1f", stage, - b.config_id, - presolve_budget_policy_name(b.policy), f.n_vars, f.n_cons, f.nnz, @@ -405,9 +163,6 @@ inline void log_presolve_budget(const char* stage, b.papilo_max_badgesize, b.probing_work_limit, b.probing_step_size, - b.intended_probe_fraction, - b.probing_work_ceiling, - (int)b.probing_ceiling_binding, b.probing_wall_limit); } diff --git a/cpp/src/mip_heuristics/presolve/probing_cache.cu b/cpp/src/mip_heuristics/presolve/probing_cache.cu index d3b7052600..546eb9f95f 100644 --- a/cpp/src/mip_heuristics/presolve/probing_cache.cu +++ b/cpp/src/mip_heuristics/presolve/probing_cache.cu @@ -890,12 +890,9 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, std::atomic problem_is_infeasible(false); size_t last_it_implied_singletons = 0; bool early_exit = false; - // Two additive terms with independently tunable weights, so a sweep can shift the budget between - // counting propagation iterations and counting probes. Both are exact counts, which is what keeps - // the budget reproducible. - const double iter_cost = (double)bound_presolve.context.settings.heuristic_params.probe_iter_work; - const double probe_cost = - (double)bound_presolve.context.settings.heuristic_params.probe_host_overhead_work; + // Two additive terms, both exact counts, which is what keeps the budget reproducible. + const double iter_cost = probing_iter_work; + const double probe_cost = probing_probe_work; // Only for the diagnostic below: a work budget buys wildly different amounts of time per // instance, so the realised rate has to be recorded to translate a budget back into seconds // afterwards. diff --git a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp index 614a23eb75..49356ff56d 100644 --- a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp +++ b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp @@ -860,7 +860,7 @@ third_party_presolve_status_t third_party_presolve_t::apply_papilo( set_presolve_parameters( papilo_presolver, category, original_n_cons, original_n_vars, max_badgesize); papilo_presolver.setVerbosityLevel(papilo::VerbosityLevel::kQuiet); - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_PAPILO_BUDGET rounds=%d badge_cap=%d tlim=%g", max_rounds, max_badgesize, time_limit); const auto papilo_t0 = std::chrono::steady_clock::now(); @@ -873,7 +873,7 @@ third_party_presolve_status_t third_party_presolve_t::apply_papilo( if (max_badgesize > 0) { effective_badge = std::min(effective_badge, max_badgesize); } // hit_tlim distinguishes "presolve converged" from "presolve was cut off mid-round", which // changes how the reduced problem below should be read. - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_PAPILO wall=%.3f tlim=%g hit_tlim=%d rounds_cap=%d badge_cap=%d badge_effective=%d", papilo_wall, time_limit, @@ -1604,16 +1604,42 @@ void papilo_postsolve_deleter::operator()(papilo::PostsolveStorage* pt delete ptr; } +template +presolve_features_t papilo_presolve_features(optimization_problem_t const& op_problem) +{ + presolve_features_t f{}; + f.n_vars = op_problem.get_n_variables(); + f.n_cons = op_problem.get_n_constraints(); + f.nnz = op_problem.get_nnz(); + + const auto var_types = op_problem.get_variable_types_host(); + const auto lower = op_problem.get_variable_lower_bounds_host(); + const auto upper = op_problem.get_variable_upper_bounds_host(); + for (size_t j = 0; j < var_types.size(); ++j) { + if (var_types[j] != var_t::INTEGER) { continue; } + f.n_int += 1.0; + if (lower[j] >= 0.0 && upper[j] <= 1.0) { f.n_bin += 1.0; } + } + + const auto offsets = op_problem.get_constraint_matrix_offsets_host(); + for (size_t i = 0; i + 1 < offsets.size(); ++i) { + f.max_row_len = std::max(f.max_row_len, offsets[i + 1] - offsets[i]); + } + return f; +} + #if MIP_INSTANTIATE_FLOAT || PDLP_INSTANTIATE_FLOAT template struct papilo_postsolve_deleter; template class third_party_presolve_t; template void papilo_round_trip(simplex::user_problem_t&); +template presolve_features_t papilo_presolve_features(optimization_problem_t const&); #endif #if MIP_INSTANTIATE_DOUBLE template struct papilo_postsolve_deleter; template class third_party_presolve_t; template void papilo_round_trip(simplex::user_problem_t&); +template presolve_features_t papilo_presolve_features(optimization_problem_t const&); #endif } // namespace cuopt::mathematical_optimization::mip diff --git a/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp b/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp index f89812e4eb..4caaae044e 100644 --- a/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp +++ b/cpp/src/mip_heuristics/presolve/third_party_presolve.hpp @@ -15,6 +15,7 @@ #include #include #include +#include #include @@ -42,6 +43,11 @@ enum class third_party_presolve_status_t { UNCHANGED, }; +// Features of the problem as the user handed it in, i.e. before any reduction. This is what Papilo +// itself will work on, so its budget is derived from these rather than from the reduced problem. +template +presolve_features_t papilo_presolve_features(optimization_problem_t const& op_problem); + template struct third_party_presolve_result_t { third_party_presolve_status_t status; diff --git a/cpp/src/mip_heuristics/solve.cu b/cpp/src/mip_heuristics/solve.cu index 1511ce1c0e..fc7b59eb68 100644 --- a/cpp/src/mip_heuristics/solve.cu +++ b/cpp/src/mip_heuristics/solve.cu @@ -333,36 +333,6 @@ mip_solution_t run_mip_solver( } } -namespace mip { - -// Features of the problem as the user handed it in, i.e. before any reduction. This is what Papilo -// itself will work on, so its budget is derived from these rather than from the reduced problem. -template -presolve_features_t papilo_presolve_features(optimization_problem_t const& op_problem) -{ - presolve_features_t f{}; - f.n_vars = op_problem.get_n_variables(); - f.n_cons = op_problem.get_n_constraints(); - f.nnz = op_problem.get_nnz(); - - const auto var_types = op_problem.get_variable_types_host(); - const auto lower = op_problem.get_variable_lower_bounds_host(); - const auto upper = op_problem.get_variable_upper_bounds_host(); - for (size_t j = 0; j < var_types.size(); ++j) { - if (var_types[j] != var_t::INTEGER) { continue; } - f.n_int += 1.0; - if (lower[j] >= 0.0 && upper[j] <= 1.0) { f.n_bin += 1.0; } - } - - const auto offsets = op_problem.get_constraint_matrix_offsets_host(); - for (size_t i = 0; i + 1 < offsets.size(); ++i) { - f.max_row_len = std::max(f.max_row_len, offsets[i + 1] - offsets[i]); - } - return f; -} - -} // namespace mip - template mip_solution_t solve_mip_helper(optimization_problem_t& op_problem, mip_solver_settings_t const& settings_const) @@ -652,7 +622,7 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p CUOPT_LOG_INFO("Papilo presolve time: %.2f", presolve_time); // What the round cap actually bought, logged here rather than inferred from the probing stage // so it is still recorded when the run never gets that far. - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_PAPILO_REDUCED nvars=%d ncons=%d nnz=%d nint=%d nbin=%d from_nvars=%.0f " "from_ncons=%.0f from_nnz=%.0f", problem.n_variables, From 1bdb2b41d96b2690e903f60eba891aa3fab2f1e0 Mon Sep 17 00:00:00 2001 From: akif Date: Tue, 4 Aug 2026 11:39:26 +0200 Subject: [PATCH 16/20] Let the work ceiling be the only bound on presolve presolve_time_ratio and presolve_max_time gave Papilo and cuOpt presolve a wall budget of min(0.1 * time_limit, 60s) each. Neither was binding at the 600s limit the set is measured at: across two repeats no Papilo call hit its tlim, the slowest finishing at 46.9s of 60 (ns1760995, whose old 243-403s blowup the badge clamp already fixed), and probing's worst was 44.8s (nw04). Both settings are removed outright, along with the 120s max_time_on_probing backstop and the dead time_ratio_of_probing_cache beside it, so the badge clamp and the probing work ceiling are what shape presolve cost while remaining_time only stops it reaching past the end of the solve. The caps did bind on short solves, where the ratio rather than the 60s term set them, so that is where this changes behaviour. Removing the settings takes their C constants, registrations, proto fields 37 and 38, and test assertions with them; the registry notes those numbers as reserved since an older client still sends them there. Every other field number is unchanged, so the wire stays compatible. Also drops the probing_wall_limit field and its log column, and the deterministic-mode branch that disabled probing when the work limit was infinite, which the previous commit made unreachable. probing_work_time_scale becomes probing_work_scale, since with no wall left it is a work coefficient and nothing converts it to seconds. CUOPT_CONFIG_ID returns with two arms so one build covers the sweep: unset or 0 is the shipping ceiling, 1 raises it 25% to measure what the extra probing buys now that no wall clips the cost of a mispredicted instance. --- .../mathematical_optimization/constants.h | 2 - .../mip/heuristics_hyper_params.hpp | 6 +- cpp/src/grpc/codegen/field_registry.yaml | 8 +-- .../codegen/generated/cuopt_remote_data.proto | 2 - .../generated_mip_settings_to_proto.inc | 2 - .../generated_proto_to_mip_settings.inc | 6 -- cpp/src/math_optimization/solver_settings.cu | 2 - .../diversity/diversity_config.hpp | 9 --- .../diversity/diversity_manager.cu | 30 +++------- .../presolve/presolve_budget_policy.hpp | 60 ++++++++++++++----- cpp/src/mip_heuristics/solve.cu | 21 +++---- cpp/src/mip_heuristics/solver.cu | 17 +++--- .../grpc/grpc_client_test.cpp | 12 ++-- cpp/tests/mip/heuristics_hyper_params_test.cu | 6 +- 14 files changed, 77 insertions(+), 106 deletions(-) diff --git a/cpp/include/cuopt/mathematical_optimization/constants.h b/cpp/include/cuopt/mathematical_optimization/constants.h index 00ef468946..04ac3355ab 100644 --- a/cpp/include/cuopt/mathematical_optimization/constants.h +++ b/cpp/include/cuopt/mathematical_optimization/constants.h @@ -98,8 +98,6 @@ #define CUOPT_MIP_HYPER_HEURISTIC_POPULATION_SIZE "mip_hyper_heuristic_population_size" #define CUOPT_MIP_HYPER_HEURISTIC_NUM_CPUFJ_THREADS "mip_hyper_heuristic_num_cpufj_threads" -#define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO "mip_hyper_heuristic_presolve_time_ratio" -#define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME "mip_hyper_heuristic_presolve_max_time" #define CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_ROUNDS "mip_hyper_heuristic_presolve_max_rounds" #define CUOPT_MIP_HYPER_HEURISTIC_PAPILO_PROBING_MAX_BADGESIZE \ "mip_hyper_heuristic_papilo_probing_max_badgesize" diff --git a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp index 94ed27f492..af72805e40 100644 --- a/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp +++ b/cpp/include/cuopt/mathematical_optimization/mip/heuristics_hyper_params.hpp @@ -19,10 +19,8 @@ namespace cuopt::mathematical_optimization { */ template struct mip_heuristics_hyper_params_t { - i_t population_size = 32; // max solutions in pool - i_t num_cpufj_threads = 8; // parallel CPU FJ climbers - f_t presolve_time_ratio = 0.1; // fraction of total time for presolve - f_t presolve_max_time = 60.0; // hard cap on presolve seconds + i_t population_size = 32; // max solutions in pool + i_t num_cpufj_threads = 8; // parallel CPU FJ climbers // Presolve budgeting. Both are derived from the problem's dimensions and structure by default // (see presolve_budget_policy.hpp); a negative value asks for that rule, and any other value diff --git a/cpp/src/grpc/codegen/field_registry.yaml b/cpp/src/grpc/codegen/field_registry.yaml index 64dcdd20e5..c13bf0db72 100644 --- a/cpp/src/grpc/codegen/field_registry.yaml +++ b/cpp/src/grpc/codegen/field_registry.yaml @@ -755,12 +755,8 @@ mip_settings: field_num: 36 type: int32 optional: true - - presolve_time_ratio: - field_num: 37 - optional: true - - presolve_max_time: - field_num: 38 - optional: true + # 37 and 38 held presolve_time_ratio and presolve_max_time, removed when presolve stopped + # taking a wall budget. Do not reuse: an older client still sends them on those numbers. - root_lp_time_ratio: field_num: 39 optional: true diff --git a/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto b/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto index 89804e2fa3..0229a25040 100644 --- a/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto +++ b/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto @@ -234,8 +234,6 @@ message MIPSolverSettings { optional double semi_continuous_big_m = 34; optional int32 population_size = 35; optional int32 num_cpufj_threads = 36; - optional double presolve_time_ratio = 37; - optional double presolve_max_time = 38; optional double root_lp_time_ratio = 39; optional double root_lp_max_time = 40; optional double rins_time_limit = 41; diff --git a/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc b/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc index 12d95aa6df..0d9919dc16 100644 --- a/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc +++ b/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc @@ -43,8 +43,6 @@ pb_settings->set_semi_continuous_big_m(settings.semi_continuous_big_m); pb_settings->set_population_size(settings.heuristic_params.population_size); pb_settings->set_num_cpufj_threads(settings.heuristic_params.num_cpufj_threads); - pb_settings->set_presolve_time_ratio(settings.heuristic_params.presolve_time_ratio); - pb_settings->set_presolve_max_time(settings.heuristic_params.presolve_max_time); pb_settings->set_root_lp_time_ratio(settings.heuristic_params.root_lp_time_ratio); pb_settings->set_root_lp_max_time(settings.heuristic_params.root_lp_max_time); pb_settings->set_rins_time_limit(settings.heuristic_params.rins_time_limit); diff --git a/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc b/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc index e14ffa1ae0..54c7238f2a 100644 --- a/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc +++ b/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc @@ -107,12 +107,6 @@ if (pb_settings.has_num_cpufj_threads()) { settings.heuristic_params.num_cpufj_threads = pb_settings.num_cpufj_threads(); } - if (pb_settings.has_presolve_time_ratio()) { - settings.heuristic_params.presolve_time_ratio = pb_settings.presolve_time_ratio(); - } - if (pb_settings.has_presolve_max_time()) { - settings.heuristic_params.presolve_max_time = pb_settings.presolve_max_time(); - } if (pb_settings.has_root_lp_time_ratio()) { settings.heuristic_params.root_lp_time_ratio = pb_settings.root_lp_time_ratio(); } diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index 181949b2d7..fc7ecc5d09 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -104,8 +104,6 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings {CUOPT_MIP_CUT_MIN_ORTHOGONALITY, &mip_settings.cut_min_orthogonality, f_t(0.0), f_t(1.0), f_t(0.5)}, {CUOPT_BARRIER_STEP_SCALE, &pdlp_settings.barrier_step_scale, f_t(0.5), f_t(0.9999), f_t(0.9)}, // MIP heuristic hyper-parameters (hidden from default --help: name contains "hyper_") - {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_TIME_RATIO, &mip_settings.heuristic_params.presolve_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for presolve"}, - {CUOPT_MIP_HYPER_HEURISTIC_PRESOLVE_MAX_TIME, &mip_settings.heuristic_params.presolve_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(60.0), "hard cap on presolve seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_TIME_RATIO, &mip_settings.heuristic_params.root_lp_time_ratio, f_t(0.0), f_t(1.0), f_t(0.1), "fraction of total time for root LP"}, {CUOPT_MIP_HYPER_HEURISTIC_ROOT_LP_MAX_TIME, &mip_settings.heuristic_params.root_lp_max_time, f_t(0.0), std::numeric_limits::infinity(), f_t(15.0), "hard cap on root LP seconds"}, {CUOPT_MIP_HYPER_HEURISTIC_RINS_TIME_LIMIT, &mip_settings.heuristic_params.rins_time_limit, f_t(0.0), std::numeric_limits::infinity(), f_t(3.0), "per-call RINS sub-MIP time"}, diff --git a/cpp/src/mip_heuristics/diversity/diversity_config.hpp b/cpp/src/mip_heuristics/diversity/diversity_config.hpp index 92d14bba40..3f71f38bc9 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_config.hpp +++ b/cpp/src/mip_heuristics/diversity/diversity_config.hpp @@ -12,15 +12,6 @@ namespace cuopt::mathematical_optimization::mip { struct diversity_config_t { - double time_ratio_of_probing_cache = 0.1; - // Backstop only. probing_work_time_scale is what bounds probing -- across 240 instances it - // stopped every run before this could fire, worst case 42.8s -- so this exists for an instance - // whose cost the proxy misses badly, not for routine use. Kept loose deliberately: tightening it - // would start shaping the common case again, which is what the work ceiling replaced. - // - // At the default presolve share (a tenth of the solve) this cannot bind below a 1200s limit, - // since the share itself is the smaller bound before then. - double max_time_on_probing = 120.0; int max_var_diff = 256; double default_time_limit = 10.; int initial_island_size = 3; diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 64699a3102..0d4f543365 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -298,18 +298,9 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ ls.constraint_prop.bounds_update.set_updated_bounds(*problem_ptr); } const auto& hp = context.settings.heuristic_params; - const bool deterministic = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC; const auto probing_features = probing_presolve_features(*problem_ptr); - auto probing_budget = evaluate_presolve_budget(hp, probing_features); + const auto probing_budget = evaluate_presolve_budget(hp, probing_features); bool run_probing_cache = !fj_only_run; - // Under the legacy policy the probing cache carries no work budget, so in deterministic mode - // nothing would bound it: the wall clock is infinite there. Keep it off, which is also what makes - // the legacy policy an unchanged baseline. - if (deterministic && - probing_budget.probing_work_limit == std::numeric_limits::infinity()) { - CUOPT_LOG_INFO("Probing-cache step disabled: deterministic mode with no work budget"); - run_probing_cache = false; - } // Allow the user to disable the probing-cache step of cuOpt's internal presolve // independently of the higher-level presolver setting. if (!context.settings.probing) { @@ -317,20 +308,13 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ run_probing_cache = false; } if (run_probing_cache) { - // Run probing cache before trivial presolve to discover variable implications. The work budget - // is what bounds this: at probing_work_time_scale it stopped every one of 240 instances before - // max_time_on_probing could fire, so the wall is a backstop for a badly mispredicted instance - // rather than the working limit it used to be. - // - // time_limit is the presolve share of the solve (presolve_time_ratio, so a tenth by default) - // and has to stay in the minimum: without it probing would be bounded only by its own wall cap, - // so a short solve would spend all of itself here -- at --time-limit 10 probing could take the - // whole 10s where it should get 1s. remaining_time only stops it reaching past the end of the - // solve. - probing_budget.probing_wall_limit = diversity_config.max_time_on_probing; + // Run probing cache before trivial presolve to discover variable implications. The work ceiling + // in presolve_budget_policy.hpp is now the only bound on its cost -- there is no wall cap and + // no presolve share -- so a proxy miss shows up as wall time rather than being clipped. Across + // 240 instances the ceiling stopped every run at a worst case of 44.7s, which is what makes + // that survivable. The timers below only stop probing reaching past the end of the solve. log_presolve_budget("PROBING", probing_features, probing_budget); - f_t time_for_probing_cache = std::min( - {(f_t)probing_budget.probing_wall_limit, time_limit, (f_t)global_timer.remaining_time()}); + f_t time_for_probing_cache = std::min(time_limit, (f_t)global_timer.remaining_time()); timer_t probing_timer{time_for_probing_cache}; const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index e8fea9565f..a59458f543 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -12,20 +12,25 @@ #include #include +#include +#include #include +#include namespace cuopt::mathematical_optimization::mip { // Work the probing loop charges per unit of effort. These are exact counts rather than timings, -// which is what makes the budget reproducible, and probing_work_time_scale below is calibrated -// against them -- changing either invalidates it. +// which is what makes the budget reproducible, and probing_work_scale below is expressed against +// them -- changing either invalidates it. inline constexpr double probing_probe_work = 0.02; // per probed variable, host overhead inline constexpr double probing_iter_work = 0.01; // per multi-probe propagation iteration -// Numerator of the probing work ceiling, divided by the cost proxy. This is the bound on probing, -// not a backstop: measured over 240 instances it stopped every run before the wall cap could fire, -// with a worst case of 44.7s, while looser scales of 4e8 and 1e9 needed the wall on 2 and 8 -// instances and spent 2-3x the total probing time to do it. +// Probing work allowed per unit of the cost proxy below; dividing by that proxy is what turns it +// into a per-instance work ceiling. It is a work coefficient, and nothing here converts it to +// seconds. Since the wall cap was removed this is the only bound on probing, which is what the +// value was picked to survive: over 240 instances it stopped every run before the 120s wall could +// fire, worst case 44.7s, while looser scales of 4e8 and 1e9 needed that wall on 2 and 8 instances +// and spent 2-3x the total probing time to do it. // // Tight enough to bound time is also tight enough to truncate, and that trade is deliberate: // probing takes its time from branch and bound, and the truncation measured neutral on solution @@ -35,7 +40,38 @@ inline constexpr double probing_iter_work = 0.01; // per multi-probe propagati // ~340x and one instance (nw04) pins the scale under every reshaping tried, including refitting the // exponent to the measured nnz^0.65. Beyond that the residual is not explained by any structural // feature; it needs throughput measured during probing rather than predicted from the problem. -inline constexpr double probing_work_time_scale = 1.5e8; +inline constexpr double probing_work_scale = 1.5e8; + +// Benchmark arm selected by CUOPT_CONFIG_ID so one build covers both points of the sweep. Unset or +// 0 is the shipping ceiling; 1 raises it 25%, which measures what the extra probing buys now that +// no wall clips the cost of a proxy miss. Read once -- the environment cannot change mid-run. +inline constexpr int n_presolve_configs = 2; +inline constexpr double probing_work_scale_arm = 1.25; + +inline double effective_probing_work_scale() +{ + static const double selected = []() -> double { + const char* raw = std::getenv("CUOPT_CONFIG_ID"); + if (raw == nullptr) { return probing_work_scale; } + int id = -1; + try { + id = std::stoi(raw); + } catch (const std::exception& e) { + CUOPT_LOG_WARN("Failed to parse CUOPT_CONFIG_ID: %s", e.what()); + return probing_work_scale; + } + if (id < 0 || id >= n_presolve_configs) { + CUOPT_LOG_WARN("CUOPT_CONFIG_ID=%d is outside [0, %d); ignoring it for presolve budgets", + id, + n_presolve_configs); + return probing_work_scale; + } + const double scale = id == 1 ? probing_work_scale * probing_work_scale_arm : probing_work_scale; + CUOPT_LOG_INFO("Using presolve budget config %d: probing work scale %.4g", id, scale); + return scale; + }(); + return selected; +} // Probed variables between work-budget checks, i.e. the granularity at which the budget can be // enforced. Work is only folded in at the step barrier, so too large a step runs unbudgeted. @@ -87,9 +123,6 @@ struct presolve_budget_t { // estimate. double probing_work_limit{std::numeric_limits::infinity()}; int probing_step_size{probing_budget_step_size}; - // Dedicated wall cap on probing. The presolve share of the solve and the global timer still apply - // on top of this; they are correctness bounds rather than tuning knobs. - double probing_wall_limit{std::numeric_limits::infinity()}; }; // Derives both presolve stages' budgets from the problem's dimensions and structure. Rounds and @@ -132,7 +165,7 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t solve_mip_helper(optimization_problem_t& op_p auto constexpr const dual_postsolve = false; if (run_presolve) { sort_csr(op_problem); - // allocate not more than 10% of the time limit to presolve. - // Note that this is not the presolve time, but the time limit for presolve. const auto& hp = settings.heuristic_params; const auto papilo_features = mip::papilo_presolve_features(op_problem); const auto papilo_budget = mip::evaluate_presolve_budget(hp, papilo_features); mip::log_presolve_budget("PAPILO", papilo_features, papilo_budget); - // The round and badge caps shape presolve, but they do not bound its cost: ns1760995 - // converges in well under 30 rounds and still spent 243-403s of a 600s budget in Papilo, - // which left no time to find a dual bound. The wall cap is what makes that survivable, so it - // stays a hard ceiling and remaining_time only stops presolve reaching past the end of the - // solve. - double presolve_time_limit = - std::min(std::min(hp.presolve_time_ratio * time_limit, hp.presolve_max_time), - timer.remaining_time()); - if (settings.determinism_mode == CUOPT_MODE_DETERMINISTIC) { - presolve_time_limit = std::numeric_limits::infinity(); - } + // Papilo carries no wall budget of its own: the badge cap is the only thing shaping its cost, + // and remaining_time just stops it reaching past the end of the solve. The removed ceiling + // was load-bearing once -- ns1760995 converges in well under 30 rounds and still spent + // 243-403s of a 600s budget here, leaving no time to find a dual bound -- so watch that + // instance and the inf-gap count if this regresses. + const double presolve_time_limit = settings.determinism_mode == CUOPT_MODE_DETERMINISTIC + ? std::numeric_limits::infinity() + : timer.remaining_time(); presolver = std::make_unique>(); auto result = presolver->apply_presolve_from_op_problem( diff --git a/cpp/src/mip_heuristics/solver.cu b/cpp/src/mip_heuristics/solver.cu index b045190e9f..f8eac0c4d8 100644 --- a/cpp/src/mip_heuristics/solver.cu +++ b/cpp/src/mip_heuristics/solver.cu @@ -207,16 +207,13 @@ solution_t mip_solver_t::run_solver() return sol; } - dm.timer = timer_; - const bool run_presolve = context.settings.presolver != presolver_t::None; - f_t time_limit = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC - ? std::numeric_limits::infinity() - : timer_.remaining_time(); - const auto& hp = context.settings.heuristic_params; - double presolve_time_limit = std::min(hp.presolve_time_ratio * time_limit, hp.presolve_max_time); - presolve_time_limit = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC - ? std::numeric_limits::infinity() - : presolve_time_limit; + dm.timer = timer_; + const bool run_presolve = context.settings.presolver != presolver_t::None; + // cuOpt presolve no longer gets a share of the solve. Its cost is bounded by the probing work + // ceiling in presolve_budget_policy.hpp, so this is only the end of the solve. + const f_t presolve_time_limit = context.settings.determinism_mode == CUOPT_MODE_DETERMINISTIC + ? std::numeric_limits::infinity() + : timer_.remaining_time(); if (std::isfinite(presolve_time_limit)) CUOPT_LOG_DEBUG("Presolve time limit: %g", presolve_time_limit); bool presolve_success = run_presolve ? dm.run_presolve(presolve_time_limit, timer_) : true; diff --git a/cpp/tests/linear_programming/grpc/grpc_client_test.cpp b/cpp/tests/linear_programming/grpc/grpc_client_test.cpp index b892e5d7ae..1e0be9d9d4 100644 --- a/cpp/tests/linear_programming/grpc/grpc_client_test.cpp +++ b/cpp/tests/linear_programming/grpc/grpc_client_test.cpp @@ -1891,8 +1891,6 @@ TEST(MapperRoundtrip, MIPSettingsAllFields) // mapping line would produce a default-valued mismatch on decode. orig.heuristic_params.population_size = 64; // default 32 orig.heuristic_params.num_cpufj_threads = 4; // default 8 - orig.heuristic_params.presolve_time_ratio = 0.2; // default 0.1 - orig.heuristic_params.presolve_max_time = 45.0; // default 60.0 orig.heuristic_params.root_lp_time_ratio = 0.25; // default 0.1 orig.heuristic_params.root_lp_max_time = 7.5; // default 15.0 orig.heuristic_params.rins_time_limit = 4.0; // default 3.0 @@ -1961,8 +1959,6 @@ TEST(MapperRoundtrip, MIPSettingsAllFields) // Heuristic hyper-parameters EXPECT_EQ(restored.heuristic_params.population_size, 64); EXPECT_EQ(restored.heuristic_params.num_cpufj_threads, 4); - EXPECT_DOUBLE_EQ(restored.heuristic_params.presolve_time_ratio, 0.2); - EXPECT_DOUBLE_EQ(restored.heuristic_params.presolve_max_time, 45.0); EXPECT_DOUBLE_EQ(restored.heuristic_params.root_lp_time_ratio, 0.25); EXPECT_DOUBLE_EQ(restored.heuristic_params.root_lp_max_time, 7.5); EXPECT_DOUBLE_EQ(restored.heuristic_params.rins_time_limit, 4.0); @@ -2496,10 +2492,10 @@ TEST(MapperRoundtrip, MIPSettingsDefaultProtoPreservesAllCppDefaults) // heuristic_params: spot-check one of each kind (int, double). EXPECT_EQ(after.heuristic_params.population_size, fresh.heuristic_params.population_size); EXPECT_EQ(after.heuristic_params.num_cpufj_threads, fresh.heuristic_params.num_cpufj_threads); - EXPECT_DOUBLE_EQ(after.heuristic_params.presolve_time_ratio, - fresh.heuristic_params.presolve_time_ratio); - EXPECT_DOUBLE_EQ(after.heuristic_params.presolve_max_time, - fresh.heuristic_params.presolve_max_time); + EXPECT_DOUBLE_EQ(after.heuristic_params.root_lp_time_ratio, + fresh.heuristic_params.root_lp_time_ratio); + EXPECT_DOUBLE_EQ(after.heuristic_params.root_lp_max_time, + fresh.heuristic_params.root_lp_max_time); EXPECT_DOUBLE_EQ(after.heuristic_params.rins_fix_rate, fresh.heuristic_params.rins_fix_rate); EXPECT_EQ(after.heuristic_params.enabled_recombiners, fresh.heuristic_params.enabled_recombiners); EXPECT_DOUBLE_EQ(after.heuristic_params.initial_infeasibility_weight, diff --git a/cpp/tests/mip/heuristics_hyper_params_test.cu b/cpp/tests/mip/heuristics_hyper_params_test.cu index c3ec4371f0..f02a5c6154 100644 --- a/cpp/tests/mip/heuristics_hyper_params_test.cu +++ b/cpp/tests/mip/heuristics_hyper_params_test.cu @@ -59,8 +59,6 @@ TEST_F(HeuristicsHyperParamsTest, CustomValuesRoundTrip) std::ofstream f(tmp_path); f << "mip_hyper_heuristic_population_size = 64\n"; f << "mip_hyper_heuristic_num_cpufj_threads = 4\n"; - f << "mip_hyper_heuristic_presolve_time_ratio = 0.2\n"; - f << "mip_hyper_heuristic_presolve_max_time = 120\n"; f << "mip_hyper_heuristic_root_lp_time_ratio = 0.05\n"; f << "mip_hyper_heuristic_root_lp_max_time = 30\n"; f << "mip_hyper_heuristic_rins_time_limit = 5\n"; @@ -82,8 +80,6 @@ TEST_F(HeuristicsHyperParamsTest, CustomValuesRoundTrip) EXPECT_EQ(hp.population_size, 64); EXPECT_EQ(hp.num_cpufj_threads, 4); - EXPECT_DOUBLE_EQ(hp.presolve_time_ratio, 0.2); - EXPECT_DOUBLE_EQ(hp.presolve_max_time, 120.0); EXPECT_DOUBLE_EQ(hp.root_lp_time_ratio, 0.05); EXPECT_DOUBLE_EQ(hp.root_lp_max_time, 30.0); EXPECT_DOUBLE_EQ(hp.rins_time_limit, 5.0); @@ -116,7 +112,7 @@ TEST_F(HeuristicsHyperParamsTest, PartialConfigKeepsDefaults) mip_heuristics_hyper_params_t defaults; EXPECT_EQ(hp.num_cpufj_threads, defaults.num_cpufj_threads); - EXPECT_DOUBLE_EQ(hp.presolve_time_ratio, defaults.presolve_time_ratio); + EXPECT_DOUBLE_EQ(hp.root_lp_time_ratio, defaults.root_lp_time_ratio); EXPECT_EQ(hp.n_of_minimums_for_exit, defaults.n_of_minimums_for_exit); EXPECT_EQ(hp.enabled_recombiners, defaults.enabled_recombiners); } From ee3465927f0de7615d25962ddf397fddd5b7d760 Mon Sep 17 00:00:00 2001 From: akif Date: Tue, 4 Aug 2026 14:59:00 +0200 Subject: [PATCH 17/20] Drop the +25% probing arm now that it measured neutral Two repeats of CUOPT_CONFIG_ID=1 against the 0b77cd173 baseline moved the set mean from 12.34 to 12.19, inside the 0.73 spread between the arm's own repeats. The arm reached only 67 of 240 instances, since work folds in at the step-128 barrier and a budget a few steps deep cannot resolve a 25% change; on those 67 it spent 23% more probing time for a median error delta of exactly 0.00, with the -0.83 mean coming from large swings that cancel in sign. The selector goes and probing_work_limit reads probing_work_scale directly again, which restores the values the default already produced. The measurement is folded into the comment on the scale, since it settles by experiment what the proxy fit had only suggested: a larger scale is not how to buy coverage back. --- .../presolve/presolve_budget_policy.hpp | 47 ++++--------------- 1 file changed, 8 insertions(+), 39 deletions(-) diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index a59458f543..7d4cea8235 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -12,10 +12,7 @@ #include #include -#include -#include #include -#include namespace cuopt::mathematical_optimization::mip { @@ -36,43 +33,15 @@ inline constexpr double probing_iter_work = 0.01; // per multi-probe propagati // probing takes its time from branch and bound, and the truncation measured neutral on solution // quality (11.78 mean error against 12.01, inside the 0.53 run-to-run noise). // -// A larger scale is not the way to buy coverage back. The proxy predicts throughput only to within -// ~340x and one instance (nw04) pins the scale under every reshaping tried, including refitting the -// exponent to the measured nnz^0.65. Beyond that the residual is not explained by any structural -// feature; it needs throughput measured during probing rather than predicted from the problem. +// A larger scale is not the way to buy coverage back, and that has now been measured rather than +// inferred: raising it 25% reached only 67 of 240 instances (the rest cannot resolve the change at +// step 128 below), spent 23% more probing time on those, and left their median error delta at +// exactly 0.00. The proxy predicts throughput only to within ~340x and one instance (nw04) pins the +// scale under every reshaping tried, including refitting the exponent to the measured nnz^0.65. +// Beyond that the residual is not explained by any structural feature; it needs throughput measured +// during probing rather than predicted from the problem. inline constexpr double probing_work_scale = 1.5e8; -// Benchmark arm selected by CUOPT_CONFIG_ID so one build covers both points of the sweep. Unset or -// 0 is the shipping ceiling; 1 raises it 25%, which measures what the extra probing buys now that -// no wall clips the cost of a proxy miss. Read once -- the environment cannot change mid-run. -inline constexpr int n_presolve_configs = 2; -inline constexpr double probing_work_scale_arm = 1.25; - -inline double effective_probing_work_scale() -{ - static const double selected = []() -> double { - const char* raw = std::getenv("CUOPT_CONFIG_ID"); - if (raw == nullptr) { return probing_work_scale; } - int id = -1; - try { - id = std::stoi(raw); - } catch (const std::exception& e) { - CUOPT_LOG_WARN("Failed to parse CUOPT_CONFIG_ID: %s", e.what()); - return probing_work_scale; - } - if (id < 0 || id >= n_presolve_configs) { - CUOPT_LOG_WARN("CUOPT_CONFIG_ID=%d is outside [0, %d); ignoring it for presolve budgets", - id, - n_presolve_configs); - return probing_work_scale; - } - const double scale = id == 1 ? probing_work_scale * probing_work_scale_arm : probing_work_scale; - CUOPT_LOG_INFO("Using presolve budget config %d: probing work scale %.4g", id, scale); - return scale; - }(); - return selected; -} - // Probed variables between work-budget checks, i.e. the granularity at which the budget can be // enforced. Work is only folded in at the step barrier, so too large a step runs unbudgeted. inline constexpr int probing_budget_step_size = 128; @@ -165,7 +134,7 @@ presolve_budget_t evaluate_presolve_budget(const mip_heuristics_hyper_params_t Date: Tue, 4 Aug 2026 15:33:36 +0200 Subject: [PATCH 18/20] Drop the presolve telemetry to DEBUG and trim its comments The budget rule is settled and measured, so the per-stage telemetry no longer needs to be on by default: PRESOLVE_BUDGET, PRESOLVE_PROBING, PRESOLVE_PROBING_WALL and the probing-cache start line all move to DEBUG. The comments recording how the ceiling was arrived at go with them, since what they argued is now either settled in the comment on probing_work_scale or recoverable from the history. Demoting is not free here. Builds set CUOPT_LOG_ACTIVE_LEVEL=INFO, so CUOPT_LOG_DEBUG compiles to nothing and diversity_manager's probing_t0 became unused-but-set under -Wall -Werror; it is marked maybe_unused. The equivalent locals in probing_cache.cu are left unannotated because they are genuinely used: probing_t0 is read by probing_wall's initializer, which still runs, and probing_wall compiles clean without a marker. Also compacts the ADAT narrowing comment to the two lines that explain the check and drops the note reserving proto fields 37 and 38. --- cpp/src/barrier/sparse_matrix_kernels.cuh | 7 ++---- cpp/src/grpc/codegen/field_registry.yaml | 2 -- .../diversity/diversity_manager.cu | 9 ++------ .../mip_heuristics/presolve/multi_probe.cuh | 6 +---- .../presolve/presolve_budget_policy.hpp | 2 +- .../mip_heuristics/presolve/probing_cache.cu | 22 +++++++------------ .../mip_heuristics/presolve/probing_cache.cuh | 5 ----- .../presolve/third_party_presolve.cpp | 2 -- cpp/src/mip_heuristics/solve.cu | 5 ----- 9 files changed, 14 insertions(+), 46 deletions(-) diff --git a/cpp/src/barrier/sparse_matrix_kernels.cuh b/cpp/src/barrier/sparse_matrix_kernels.cuh index 33fc2bfe41..0ce8447307 100644 --- a/cpp/src/barrier/sparse_matrix_kernels.cuh +++ b/cpp/src/barrier/sparse_matrix_kernels.cuh @@ -145,11 +145,8 @@ void multiply_kernels(raft::handle_t const* handle, int64_t ADAT_num_rows, ADAT_num_cols, ADAT_nnz1; RAFT_CUSPARSE_TRY( cusparseSpMatGetSize(cusparse_data.matADAT_descr, &ADAT_num_rows, &ADAT_num_cols, &ADAT_nnz1)); - // cuSPARSE sizes the product in 64 bits, but the CSR arrays are indexed by i_t. A tall problem - // with a few near-full columns leaves ADAT dense enough to pass that range even after dense - // columns are eliminated, and narrowing would hand RMM a negative count that resurfaces as an - // unrelated "size overflows device_uvector storage" from deep inside the allocator. Report the - // real cause instead, as the capacity failure the caller already knows how to fall back from. + // cuSPARSE sizes the product in 64 bits while the CSR arrays are indexed by i_t; narrowing would + // reach RMM as a negative count and surface as an unrelated device_uvector overflow. if (ADAT_nnz1 > std::numeric_limits::max()) { throw rmm::out_of_memory( "ADAT needs " + std::to_string(ADAT_nnz1) + " nonzeros over " + diff --git a/cpp/src/grpc/codegen/field_registry.yaml b/cpp/src/grpc/codegen/field_registry.yaml index c13bf0db72..10fe5c0c7a 100644 --- a/cpp/src/grpc/codegen/field_registry.yaml +++ b/cpp/src/grpc/codegen/field_registry.yaml @@ -755,8 +755,6 @@ mip_settings: field_num: 36 type: int32 optional: true - # 37 and 38 held presolve_time_ratio and presolve_max_time, removed when presolve stopped - # taking a wall budget. Do not reuse: an older client still sends them on those numbers. - root_lp_time_ratio: field_num: 39 optional: true diff --git a/cpp/src/mip_heuristics/diversity/diversity_manager.cu b/cpp/src/mip_heuristics/diversity/diversity_manager.cu index 0d4f543365..9d70ae17ee 100644 --- a/cpp/src/mip_heuristics/diversity/diversity_manager.cu +++ b/cpp/src/mip_heuristics/diversity/diversity_manager.cu @@ -308,15 +308,10 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ run_probing_cache = false; } if (run_probing_cache) { - // Run probing cache before trivial presolve to discover variable implications. The work ceiling - // in presolve_budget_policy.hpp is now the only bound on its cost -- there is no wall cap and - // no presolve share -- so a proxy miss shows up as wall time rather than being clipped. Across - // 240 instances the ceiling stopped every run at a worst case of 44.7s, which is what makes - // that survivable. The timers below only stop probing reaching past the end of the solve. log_presolve_budget("PROBING", probing_features, probing_budget); f_t time_for_probing_cache = std::min(time_limit, (f_t)global_timer.remaining_time()); timer_t probing_timer{time_for_probing_cache}; - const auto probing_t0 = std::chrono::steady_clock::now(); + [[maybe_unused]] const auto probing_t0 = std::chrono::steady_clock::now(); // this function computes probing cache, finds singletons, substitutions and changes the problem bool problem_is_infeasible = compute_probing_cache(ls.constraint_prop.bounds_update, *problem_ptr, @@ -324,7 +319,7 @@ bool diversity_manager_t::run_presolve(f_t time_limit, timer_t global_ probing_budget.probing_work_limit, (size_t)probing_budget.probing_step_size); problem_ptr->handle_ptr->sync_stream(); - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_PROBING_WALL wall=%.3f", std::chrono::duration(std::chrono::steady_clock::now() - probing_t0).count()); if (problem_is_infeasible) { return false; } diff --git a/cpp/src/mip_heuristics/presolve/multi_probe.cuh b/cpp/src/mip_heuristics/presolve/multi_probe.cuh index 2f2f2f038d..0dfa780694 100644 --- a/cpp/src/mip_heuristics/presolve/multi_probe.cuh +++ b/cpp/src/mip_heuristics/presolve/multi_probe.cuh @@ -75,11 +75,7 @@ class multi_probe_t { bool skip_0; bool skip_1; settings_t settings; - // When set, the number of propagation iterations run is accumulated here rather than on any - // shared counter. The probing cache runs one multi_probe_t per OMP task, so a shared counter - // would be an unsynchronized read-modify-write and would make the resulting budget - // nondeterministic. The owner folds the per-task counts in at a barrier, in a fixed order, and - // applies its own cost model. + // Per-task iteration count; a shared counter would race and make the budget nondeterministic. double* local_iter_accumulator = nullptr; bool compute_stats = true; bool init_changed_constraints = true; diff --git a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp index 7d4cea8235..1b18fdec96 100644 --- a/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp +++ b/cpp/src/mip_heuristics/presolve/presolve_budget_policy.hpp @@ -145,7 +145,7 @@ inline void log_presolve_budget(const char* stage, const presolve_features_t& f, const presolve_budget_t& b) { - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_BUDGET stage=%s nvars=%.0f ncons=%.0f nnz=%.0f nint=%.0f " "nbin=%.0f arl=%.3f acl=%.3f maxrow=%.0f density=%.3e intfrac=%.3f binfrac=%.3f " "rounds=%d badge=%d work=%.3f step=%d", diff --git a/cpp/src/mip_heuristics/presolve/probing_cache.cu b/cpp/src/mip_heuristics/presolve/probing_cache.cu index 546eb9f95f..38d076581f 100644 --- a/cpp/src/mip_heuristics/presolve/probing_cache.cu +++ b/cpp/src/mip_heuristics/presolve/probing_cache.cu @@ -872,8 +872,6 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, std::vector>> modification_vector_pool(num_tasks); std::vector>> substitution_vector_pool(num_tasks); - // Each task counts its own propagation iterations; the counts are folded in at the step barrier - // below, where only one thread is running. std::vector iter_accum_pool(num_tasks, 0.0); // Initialize multi_probe_presolve_pool @@ -890,16 +888,12 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, std::atomic problem_is_infeasible(false); size_t last_it_implied_singletons = 0; bool early_exit = false; - // Two additive terms, both exact counts, which is what keeps the budget reproducible. - const double iter_cost = probing_iter_work; - const double probe_cost = probing_probe_work; - // Only for the diagnostic below: a work budget buys wildly different amounts of time per - // instance, so the realised rate has to be recorded to translate a budget back into seconds - // afterwards. - const auto probing_t0 = std::chrono::steady_clock::now(); - double iters_done = 0.0; - size_t probes_done = 0; - double work_used = 0.0; + const double iter_cost = probing_iter_work; + const double probe_cost = probing_probe_work; + const auto probing_t0 = std::chrono::steady_clock::now(); + double iters_done = 0.0; + size_t probes_done = 0; + double work_used = 0.0; // Work is only folded in at the step barrier, so the step size is also the granularity at which // the budget can be enforced: too large and a single step runs effectively unbudgeted. const size_t step_size = min(step_size_hint, priority_indices.size()); @@ -909,7 +903,7 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, // are visible before any per-thread kernel can reference that memory. problem.handle_ptr->sync_stream(); - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "Running probing cache with %zu tasks (%zu candidate vars, work limit %.3f, step %zu)", num_tasks, priority_indices.size(), @@ -978,7 +972,7 @@ bool compute_probing_cache(bound_presolve_t& bound_presolve, apply_substitution_queue_to_problem(substitution_vector_pool, problem); const double probing_wall = std::chrono::duration(std::chrono::steady_clock::now() - probing_t0).count(); - CUOPT_LOG_INFO( + CUOPT_LOG_DEBUG( "PRESOLVE_PROBING probes=%zu candidates=%zu iters=%.0f work=%.3f work_limit=%.3f step=%zu " "iter_cost=%.5f probe_cost=%.5f wall=%.3f wall_limit=%.3f units_per_s=%.1f " "budget_exhausted=%d early_exit=%d timed_out=%d cached=%lu implied_singletons=%lu", diff --git a/cpp/src/mip_heuristics/presolve/probing_cache.cuh b/cpp/src/mip_heuristics/presolve/probing_cache.cuh index 74526d1ba2..24d9a9cfc1 100644 --- a/cpp/src/mip_heuristics/presolve/probing_cache.cuh +++ b/cpp/src/mip_heuristics/presolve/probing_cache.cuh @@ -123,8 +123,6 @@ class lb_probing_cache_t { std::unordered_map, 2>> probing_cache; }; -// Features of the Papilo-reduced problem, which is what the probing cache actually runs on. These -// differ from the features Papilo's own budget was derived from, sometimes by a lot. template presolve_features_t probing_presolve_features(problem_t const& problem) { @@ -142,9 +140,6 @@ presolve_features_t probing_presolve_features(problem_t const& problem return f; } -// `work_limit` bounds probing in work units, checked at every step barrier and therefore -// independent of thread count and wall clock. `step_size_hint` is the number of variables probed -// per step, i.e. the granularity at which the budget can be enforced. template bool compute_probing_cache(bound_presolve_t& bound_presolve, problem_t& problem, diff --git a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp index 49356ff56d..8f0c62c053 100644 --- a/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp +++ b/cpp/src/mip_heuristics/presolve/third_party_presolve.cpp @@ -716,8 +716,6 @@ void set_presolve_options(papilo::Presolve& presolver, presolver.getPresolveOptions().tlim = time_limit; presolver.getPresolveOptions().threads = num_cpu_threads; // user setting or 0 (automatic) presolver.getPresolveOptions().feastol = 1e-5; - // A round cap bounds presolve independently of the clock, which is the only thing that bounds it - // in deterministic mode where tlim is infinite. <=0 keeps Papilo's default of unlimited rounds. if (max_rounds > 0) { presolver.getPresolveOptions().maxrounds = max_rounds; } if (dual_postsolve) { presolver.getPresolveOptions().componentsmaxint = -1; diff --git a/cpp/src/mip_heuristics/solve.cu b/cpp/src/mip_heuristics/solve.cu index 33d01485a8..f8658237fd 100644 --- a/cpp/src/mip_heuristics/solve.cu +++ b/cpp/src/mip_heuristics/solve.cu @@ -565,11 +565,6 @@ mip_solution_t solve_mip_helper(optimization_problem_t& op_p const auto papilo_budget = mip::evaluate_presolve_budget(hp, papilo_features); mip::log_presolve_budget("PAPILO", papilo_features, papilo_budget); - // Papilo carries no wall budget of its own: the badge cap is the only thing shaping its cost, - // and remaining_time just stops it reaching past the end of the solve. The removed ceiling - // was load-bearing once -- ns1760995 converges in well under 30 rounds and still spent - // 243-403s of a 600s budget here, leaving no time to find a dual bound -- so watch that - // instance and the inf-gap count if this regresses. const double presolve_time_limit = settings.determinism_mode == CUOPT_MODE_DETERMINISTIC ? std::numeric_limits::infinity() : timer.remaining_time(); From ec3adedd1e0189e691541033084d6a38a71dca63 Mon Sep 17 00:00:00 2001 From: akif Date: Tue, 4 Aug 2026 18:07:36 +0200 Subject: [PATCH 19/20] fix ai comments --- cpp/src/grpc/codegen/field_registry.yaml | 11 +++++++++++ .../grpc/codegen/generated/cuopt_remote_data.proto | 2 ++ .../generated/generated_mip_settings_to_proto.inc | 2 ++ .../generated/generated_proto_to_mip_settings.inc | 6 ++++++ .../linear_programming/grpc/grpc_client_test.cpp | 7 +++++++ 5 files changed, 28 insertions(+) diff --git a/cpp/src/grpc/codegen/field_registry.yaml b/cpp/src/grpc/codegen/field_registry.yaml index 10fe5c0c7a..2fb7896027 100644 --- a/cpp/src/grpc/codegen/field_registry.yaml +++ b/cpp/src/grpc/codegen/field_registry.yaml @@ -755,6 +755,9 @@ mip_settings: field_num: 36 type: int32 optional: true + # 37 and 38 held presolve_time_ratio and presolve_max_time, removed when + # presolve stopped taking a wall budget. Do not reuse: an older client + # still sends them on those numbers. - root_lp_time_ratio: field_num: 39 optional: true @@ -799,6 +802,14 @@ mip_settings: - related_vars_time_limit: field_num: 51 optional: true + - presolve_max_rounds: + field_num: 53 + type: int32 + optional: true + - papilo_probing_max_badgesize: + field_num: 54 + type: int32 + optional: true # ───────────────────────────────────────────────────────────────────────────── # Optimization Problem (cpu_optimization_problem_t) diff --git a/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto b/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto index 0229a25040..4b1e36d134 100644 --- a/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto +++ b/cpp/src/grpc/codegen/generated/cuopt_remote_data.proto @@ -248,6 +248,8 @@ message MIPSolverSettings { optional double relaxed_lp_time_limit = 50; optional double related_vars_time_limit = 51; optional int32 zero_half_cuts = 52; + optional int32 presolve_max_rounds = 53; + optional int32 papilo_probing_max_badgesize = 54; } message PDLPWarmStartData { diff --git a/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc b/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc index 0d9919dc16..e92e156053 100644 --- a/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc +++ b/cpp/src/grpc/codegen/generated/generated_mip_settings_to_proto.inc @@ -56,3 +56,5 @@ pb_settings->set_cycle_detection_length(settings.heuristic_params.cycle_detection_length); pb_settings->set_relaxed_lp_time_limit(settings.heuristic_params.relaxed_lp_time_limit); pb_settings->set_related_vars_time_limit(settings.heuristic_params.related_vars_time_limit); + pb_settings->set_presolve_max_rounds(settings.heuristic_params.presolve_max_rounds); + pb_settings->set_papilo_probing_max_badgesize(settings.heuristic_params.papilo_probing_max_badgesize); diff --git a/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc b/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc index 54c7238f2a..04241b3407 100644 --- a/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc +++ b/cpp/src/grpc/codegen/generated/generated_proto_to_mip_settings.inc @@ -146,3 +146,9 @@ if (pb_settings.has_related_vars_time_limit()) { settings.heuristic_params.related_vars_time_limit = pb_settings.related_vars_time_limit(); } + if (pb_settings.has_presolve_max_rounds()) { + settings.heuristic_params.presolve_max_rounds = pb_settings.presolve_max_rounds(); + } + if (pb_settings.has_papilo_probing_max_badgesize()) { + settings.heuristic_params.papilo_probing_max_badgesize = pb_settings.papilo_probing_max_badgesize(); + } diff --git a/cpp/tests/linear_programming/grpc/grpc_client_test.cpp b/cpp/tests/linear_programming/grpc/grpc_client_test.cpp index 1e0be9d9d4..f8fed6eee3 100644 --- a/cpp/tests/linear_programming/grpc/grpc_client_test.cpp +++ b/cpp/tests/linear_programming/grpc/grpc_client_test.cpp @@ -1891,6 +1891,8 @@ TEST(MapperRoundtrip, MIPSettingsAllFields) // mapping line would produce a default-valued mismatch on decode. orig.heuristic_params.population_size = 64; // default 32 orig.heuristic_params.num_cpufj_threads = 4; // default 8 + orig.heuristic_params.presolve_max_rounds = 12; // default -1 + orig.heuristic_params.papilo_probing_max_badgesize = 64; // default -1 orig.heuristic_params.root_lp_time_ratio = 0.25; // default 0.1 orig.heuristic_params.root_lp_max_time = 7.5; // default 15.0 orig.heuristic_params.rins_time_limit = 4.0; // default 3.0 @@ -1959,6 +1961,8 @@ TEST(MapperRoundtrip, MIPSettingsAllFields) // Heuristic hyper-parameters EXPECT_EQ(restored.heuristic_params.population_size, 64); EXPECT_EQ(restored.heuristic_params.num_cpufj_threads, 4); + EXPECT_EQ(restored.heuristic_params.presolve_max_rounds, 12); + EXPECT_EQ(restored.heuristic_params.papilo_probing_max_badgesize, 64); EXPECT_DOUBLE_EQ(restored.heuristic_params.root_lp_time_ratio, 0.25); EXPECT_DOUBLE_EQ(restored.heuristic_params.root_lp_max_time, 7.5); EXPECT_DOUBLE_EQ(restored.heuristic_params.rins_time_limit, 4.0); @@ -2492,6 +2496,9 @@ TEST(MapperRoundtrip, MIPSettingsDefaultProtoPreservesAllCppDefaults) // heuristic_params: spot-check one of each kind (int, double). EXPECT_EQ(after.heuristic_params.population_size, fresh.heuristic_params.population_size); EXPECT_EQ(after.heuristic_params.num_cpufj_threads, fresh.heuristic_params.num_cpufj_threads); + EXPECT_EQ(after.heuristic_params.presolve_max_rounds, fresh.heuristic_params.presolve_max_rounds); + EXPECT_EQ(after.heuristic_params.papilo_probing_max_badgesize, + fresh.heuristic_params.papilo_probing_max_badgesize); EXPECT_DOUBLE_EQ(after.heuristic_params.root_lp_time_ratio, fresh.heuristic_params.root_lp_time_ratio); EXPECT_DOUBLE_EQ(after.heuristic_params.root_lp_max_time, From 9fd5e1e7e83f58b9ded764c7f0c2679c2d197bdb Mon Sep 17 00:00:00 2001 From: akif Date: Mon, 10 Aug 2026 15:27:02 +0200 Subject: [PATCH 20/20] remove leftover file --- cpp/src/cuts/zero_half_mod2.cpp | 727 -------------------------------- 1 file changed, 727 deletions(-) delete mode 100644 cpp/src/cuts/zero_half_mod2.cpp diff --git a/cpp/src/cuts/zero_half_mod2.cpp b/cpp/src/cuts/zero_half_mod2.cpp deleted file mode 100644 index ac55131563..0000000000 --- a/cpp/src/cuts/zero_half_mod2.cpp +++ /dev/null @@ -1,727 +0,0 @@ -/* clang-format off */ -/* - * SPDX-FileCopyrightText: Copyright (c) 2026, NVIDIA CORPORATION & AFFILIATES. All rights reserved. - * SPDX-License-Identifier: Apache-2.0 - */ -/* clang-format on */ - -#include - -#include -#include - -#include -#include -#include -#include -#include -#include -#include -#include - -namespace cuopt::mathematical_optimization::mip { - -using simplex::lp_problem_t; -using simplex::simplex_solver_settings_t; -using simplex::variable_type_t; - -namespace { - -template -void symmetric_difference_sorted(const std::vector& a, - const std::vector& b, - std::vector& result) -{ - result.clear(); - result.reserve(a.size() + b.size()); - std::set_symmetric_difference(a.begin(), a.end(), b.begin(), b.end(), std::back_inserter(result)); -} - -template -struct mod2_parity_row_t { - std::vector parity; - bool rhs_parity{false}; -}; - -template -struct mod2_candidate_t : mod2_parity_row_t { - inequality_t transformed_inequality; - bool reversible{false}; -}; - -template -struct mod2_basis_row_t { - std::vector parity; - std::vector combination; - bool rhs{false}; -}; - -template -struct mod2_row_order_t { - const std::vector& rows; - - bool operator()(i_t a, i_t b) const - { - if (rows[a].parity.size() != rows[b].parity.size()) { - return rows[a].parity.size() < rows[b].parity.size(); - } - return rows[a].rhs_parity < rows[b].rhs_parity; - } -}; - -template -std::vector> find_mod2_row_combinations(const std::vector& rows, - i_t max_combination_size, - i_t max_combinations, - f_t* work_estimate, - f_t max_work_estimate) -{ - cuopt_assert(max_combination_size > 0, "Maximum GF(2) combination size must be positive"); - cuopt_assert(max_combinations > 0, "Maximum number of GF(2) combinations must be positive"); - - i_t max_index = -1; - f_t input_scan_work = 0.0; - for (const auto& row : rows) { - input_scan_work += (f_t)(row.parity.size() + 1); - cuopt_assert(std::is_sorted(row.parity.begin(), row.parity.end()), - "GF(2) parity rows must be sorted"); - cuopt_assert(std::adjacent_find(row.parity.begin(), row.parity.end()) == row.parity.end(), - "GF(2) parity rows must not contain duplicates"); - if (!row.parity.empty()) { - cuopt_assert(row.parity.front() >= 0, "GF(2) parity index must be nonnegative"); - max_index = std::max(max_index, row.parity.back()); - } - } - if (add_work_estimate(input_scan_work, work_estimate, max_work_estimate)) { return {}; } - - std::vector permutation(rows.size()); - std::iota(permutation.begin(), permutation.end(), 0); - const f_t sort_work = (f_t)permutation.size() * std::log2((f_t)permutation.size() + (f_t)1.0); - if (add_work_estimate(sort_work, work_estimate, max_work_estimate)) { return {}; } - // this is to process small/sparse rows first, for faster perf and smaller combinations - std::stable_sort(permutation.begin(), permutation.end(), mod2_row_order_t{rows}); - - if (add_work_estimate((f_t)(max_index + 1), work_estimate, max_work_estimate)) { return {}; } - std::vector pivot_to_basis((size_t)(max_index + 1), -1); - std::vector> basis; - basis.reserve(std::min(rows.size(), (size_t)(max_index + 1))); - std::vector> combinations; - combinations.reserve(std::min((size_t)max_combinations, rows.size())); - - std::vector parity_tmp; - std::vector combination_tmp; - for (const i_t candidate : permutation) { - f_t candidate_work = (f_t)(rows[candidate].parity.size() + 2); - mod2_basis_row_t current; - current.parity = rows[candidate].parity; - current.combination = {candidate}; - current.rhs = rows[candidate].rhs_parity; - - bool abandoned = false; - while (!current.parity.empty()) { - const i_t pivot = current.parity.front(); - const i_t basis_index = pivot_to_basis[pivot]; - // pivot has not been seen before - if (basis_index < 0) { break; } - - const auto& pivot_row = basis[basis_index]; - candidate_work += (f_t)(current.parity.size() + pivot_row.parity.size() + - current.combination.size() + pivot_row.combination.size()); - symmetric_difference_sorted(current.parity, pivot_row.parity, parity_tmp); - symmetric_difference_sorted(current.combination, pivot_row.combination, combination_tmp); - if (combination_tmp.size() > (size_t)max_combination_size) { - abandoned = true; - break; - } - current.parity.swap(parity_tmp); - current.combination.swap(combination_tmp); - current.rhs = current.rhs != pivot_row.rhs; - } - if (add_work_estimate(candidate_work, work_estimate, max_work_estimate)) { break; } - if (abandoned) { continue; } - - // when reduced, add to combinations and continue, don't add to basis - if (current.parity.empty()) { - if (current.rhs && !current.combination.empty()) { - combinations.push_back(std::move(current.combination)); - if (combinations.size() >= (size_t)max_combinations) { break; } - } - continue; - } - - const i_t pivot = current.parity.front(); - pivot_to_basis[pivot] = (i_t)basis.size(); - basis.push_back(std::move(current)); - } - return combinations; -} - -template -i_t mod2_integral_scale(const inequality_t& inequality, - const std::vector& var_types, - const std::vector& transformed_xstar, - i_t max_integral_scale, - f_t row_tight_tol, - f_t coefficient_integral_tol, - f_t start_time, - f_t time_limit, - f_t& work_estimate) -{ - if (toc(start_time) >= time_limit) { return i_t{0}; } - f_t scale_work = 0.0; - for (i_t scale = 1; scale <= max_integral_scale; ++scale) { - scale_work += 1.0; - bool integral = true; - const f_t scaled_rhs = (f_t)scale * inequality.rhs; - if (std::abs(scaled_rhs - std::round(scaled_rhs)) > - coefficient_integral_tol * std::max((f_t)1.0, std::abs(scaled_rhs))) { - integral = false; - } - for (i_t k = 0; integral && k < (i_t)inequality.size(); ++k) { - scale_work += 1.0; - const i_t j = inequality.index(k); - if (var_types[j] == variable_type_t::CONTINUOUS || transformed_xstar[j] <= row_tight_tol) { - continue; - } - const f_t scaled_coefficient = (f_t)scale * inequality.coeff(k); - if (std::abs(scaled_coefficient - std::round(scaled_coefficient)) > - coefficient_integral_tol * std::max((f_t)1.0, std::abs(scaled_coefficient))) { - integral = false; - } - } - if (integral) { - work_estimate += scale_work; - return scale; - } - } - work_estimate += scale_work; - return i_t{0}; -} - -template -std::vector> mod2_collect_candidates( - complemented_mixed_integer_rounding_cut_t& complemented_mir, - const lp_problem_t& lp, - csr_matrix_t& Arow, - const variable_bounds_t& variable_bounds, - const std::vector& var_types, - const std::vector& transformed_xstar, - f_t start_time, - f_t time_limit, - f_t& work_estimate, - f_t max_work_estimate, - bool& work_limit_reached) -{ - constexpr i_t max_integral_scale = 1000; - const i_t max_integer_row_length = 1000 + lp.num_cols / 10; - constexpr f_t row_tight_tol = (f_t)1e-6; - constexpr f_t coefficient_integral_tol = (f_t)1e-6; - - std::vector> candidates; - candidates.reserve(lp.num_rows); - for (i_t row = 0; row < lp.num_rows; ++row) { - if (toc(start_time) >= time_limit || work_limit_reached) { break; } - const i_t slack = complemented_mir.slack_cols(row); - if (slack < 0 || transformed_xstar[slack] > row_tight_tol) { continue; } - - const i_t row_length = Arow.row_start[row + 1] - Arow.row_start[row]; - if (row_length > max_integer_row_length) { continue; } - const f_t row_work = - (f_t)(8 * row_length + 5) + (f_t)row_length * std::log2((f_t)row_length + (f_t)1.0); - if (add_work_estimate(row_work, &work_estimate, max_work_estimate, &work_limit_reached)) { - break; - } - inequality_t inequality(Arow, row, lp.rhs[row]); - complemented_mir.transform_inequality(variable_bounds, var_types, inequality); - inequality.sort(); - - // Every LP row is an equality after slack insertion. Remove a zero-valued transformed slack - // in the direction that preserves a valid >= inequality. - i_t slack_position = -1; - for (i_t k = 0; k < (i_t)inequality.size(); ++k) { - if (inequality.index(k) == slack) { - slack_position = k; - break; - } - } - if (slack_position < 0 || inequality.coeff(slack_position) == 0.0) { continue; } - // we want a row that is a.x >= b - if (inequality.coeff(slack_position) > 0.0) { inequality.negate(); } - inequality.vector.x[slack_position] = 0.0; - inequality_t squeezed_inequality(lp.num_cols); - inequality.squeeze(squeezed_inequality); - inequality = std::move(squeezed_inequality); - - // Continuous variables must be at their selected bounds to participate in the parity system. - bool continuous_at_bounds = true; - for (i_t k = 0; k < (i_t)inequality.size(); ++k) { - const i_t j = inequality.index(k); - if (var_types[j] == variable_type_t::CONTINUOUS && - std::abs(inequality.coeff(k)) > coefficient_integral_tol && - transformed_xstar[j] > row_tight_tol) { - continuous_at_bounds = false; - break; - } - } - if (!continuous_at_bounds) { continue; } - - const i_t scale = mod2_integral_scale(inequality, - var_types, - transformed_xstar, - max_integral_scale, - row_tight_tol, - coefficient_integral_tol, - start_time, - time_limit, - work_estimate); - if (work_estimate > max_work_estimate) { - work_limit_reached = true; - break; - } - // no integral scale found or time limit reached - if (scale == 0) { continue; } - if (scale != 1) { inequality.scale((f_t)scale); } - - mod2_candidate_t candidate; - candidate.transformed_inequality = std::move(inequality); - candidate.rhs_parity = (std::abs(std::llround(candidate.transformed_inequality.rhs)) % 2) != 0; - // checks if this could be safely reversed - candidate.reversible = std::abs(lp.upper[slack] - lp.lower[slack]) <= row_tight_tol; - for (i_t k = 0; k < (i_t)candidate.transformed_inequality.size(); ++k) { - const i_t j = candidate.transformed_inequality.index(k); - if (var_types[j] == variable_type_t::CONTINUOUS || transformed_xstar[j] <= row_tight_tol) { - continue; - } - const auto coefficient = std::llround(candidate.transformed_inequality.coeff(k)); - if ((std::abs(coefficient) % 2) != 0) { candidate.parity.push_back(j); } - } - if (candidate.parity.size() > (size_t)max_integer_row_length) { continue; } - candidates.push_back(std::move(candidate)); - } - return candidates; -} - -template -void mod2_add_transformed_zero_half_cut( - complemented_mixed_integer_rounding_cut_t& complemented_mir, - cut_pool_t& cut_pool, - const lp_problem_t& lp, - csr_matrix_t& Arow, - const variable_bounds_t& variable_bounds, - const std::vector& var_types, - const std::vector& xstar, - inequality_t transformed_cut, - f_t min_violation, - f_t& work_estimate, - i_t& cuts_added) -{ - work_estimate += (f_t)(10 * transformed_cut.size() + 1); - complemented_mir.untransform_inequality(variable_bounds, var_types, transformed_cut); - complemented_mir.remove_small_coefficients(lp.lower, lp.upper, transformed_cut); - complemented_mir.substitute_slacks(lp, Arow, transformed_cut); - complemented_mir.remove_small_coefficients(lp.lower, lp.upper, transformed_cut); - const f_t violation = complemented_mir.compute_violation(transformed_cut, xstar); - if (violation > min_violation) { - cut_pool.add_cut(cut_type_t::ZERO_HALF, transformed_cut); - ++cuts_added; - } -} - -template -void mod2_generate_cuts_from_aggregate( - complemented_mixed_integer_rounding_cut_t& complemented_mir, - cut_pool_t& cut_pool, - const lp_problem_t& lp, - const simplex_solver_settings_t& settings, - csr_matrix_t& Arow, - const variable_bounds_t& variable_bounds, - const std::vector& var_types, - const std::vector& xstar, - const std::vector& transformed_xstar, - const inequality_t& oriented_aggregate, - f_t min_violation, - f_t start_time, - f_t& work_estimate, - f_t max_work_estimate, - bool& work_limit_reached, - i_t& cuts_added) -{ - work_estimate += (f_t)(3 * oriented_aggregate.size() + 1); - inequality_t mir_cut(lp.num_cols); - const bool mir_cut_generated = complemented_mir.generate_cut_nonnegative_maintain_indicies( - oriented_aggregate, var_types, mir_cut); - if (mir_cut_generated) { - mod2_add_transformed_zero_half_cut(complemented_mir, - cut_pool, - lp, - Arow, - variable_bounds, - var_types, - xstar, - std::move(mir_cut), - min_violation, - work_estimate, - cuts_added); - } - - if (work_estimate > max_work_estimate) { - work_limit_reached = true; - return; - } - inequality_t lifted_cover_cut(lp.num_cols); - bool lifted_cover_cut_generated = false; - if (toc(start_time) < settings.time_limit) { - lifted_cover_cut_generated = - complemented_mir.generate_lifted_mixed_binary_cover(oriented_aggregate, - var_types, - transformed_xstar, - lifted_cover_cut, - work_estimate, - max_work_estimate); - } - if (lifted_cover_cut_generated) { - mod2_add_transformed_zero_half_cut(complemented_mir, - cut_pool, - lp, - Arow, - variable_bounds, - var_types, - xstar, - std::move(lifted_cover_cut), - min_violation, - work_estimate, - cuts_added); - } - if (work_estimate > max_work_estimate) { work_limit_reached = true; } -} - -template -struct lifted_cover_order_t { - const std::vector& solution_value; - const inequality_t& base; - f_t tolerance; - - bool operator()(int a, int b) const - { - const bool a_at_upper = solution_value[a] >= 1.0 - tolerance; - const bool b_at_upper = solution_value[b] >= 1.0 - tolerance; - if (a_at_upper != b_at_upper) { return a_at_upper; } - const f_t contribution_a = solution_value[a] * base.coeff(a); - const f_t contribution_b = solution_value[b] * base.coeff(b); - if (contribution_a != contribution_b) { return contribution_a > contribution_b; } - return base.coeff(a) > base.coeff(b); - } -}; - -template -f_t lifted_cover_coefficient( - f_t coefficient, const std::vector& prefix, size_t p, f_t lambda, f_t tolerance) -{ - for (size_t h = 0; h < p; ++h) { - if (coefficient <= prefix[h] - lambda + tolerance) { return (f_t)h * lambda; } - if (coefficient <= prefix[h] + tolerance) { - return (f_t)(h + 1) * lambda + coefficient - prefix[h]; - } - } - return (f_t)p * lambda + coefficient - prefix[p - 1]; -} - -} // namespace - -std::vector> find_mod2_row_combinations_for_test( - const std::vector>& parity_rows, - const std::vector& rhs_parity, - int max_combination_size, - int max_combinations) -{ - return find_mod2_row_combinations_for_test(parity_rows, - rhs_parity, - max_combination_size, - max_combinations, - std::numeric_limits::infinity(), - nullptr); -} - -std::vector> find_mod2_row_combinations_for_test( - const std::vector>& parity_rows, - const std::vector& rhs_parity, - int max_combination_size, - int max_combinations, - double max_work_estimate, - double* work_estimate_out) -{ - cuopt_assert(parity_rows.size() == rhs_parity.size(), - "GF(2) parity row and rhs sizes must match"); - std::vector> rows; - rows.reserve(parity_rows.size()); - for (size_t i = 0; i < parity_rows.size(); ++i) { - rows.push_back({parity_rows[i], rhs_parity[i] != 0}); - } - - double work_estimate = 0.0; - auto combinations = find_mod2_row_combinations( - rows, max_combination_size, max_combinations, &work_estimate, max_work_estimate); - if (work_estimate_out != nullptr) { *work_estimate_out = work_estimate; } - return combinations; -} - -template -bool generate_mod2_zero_half_cuts(cut_pool_t& cut_pool, - const lp_problem_t& lp, - const simplex_solver_settings_t& settings, - csr_matrix_t& Arow, - const std::vector& new_slacks, - const std::vector& var_types, - const std::vector& xstar, - variable_bounds_t& variable_bounds, - f_t start_time, - f_t& work_estimate) -{ - constexpr i_t max_combination_size = 64; - constexpr i_t max_row_combinations = 1000; - constexpr f_t min_violation = (f_t)1e-6; - const f_t max_work_estimate = work_estimate + (f_t)1e8; - bool work_limit_reached = false; - - if (add_work_estimate((f_t)(3 * lp.num_cols) + (f_t)(variable_bounds.upper_variables.size() + - variable_bounds.lower_variables.size()), - &work_estimate, - max_work_estimate, - &work_limit_reached)) { - return false; - } - complemented_mixed_integer_rounding_cut_t complemented_mir(lp, settings, new_slacks); - std::vector transformed_xstar; - complemented_mir.bound_substitution( - lp, variable_bounds, var_types, xstar, transformed_xstar, true); - - auto candidates = mod2_collect_candidates(complemented_mir, - lp, - Arow, - variable_bounds, - var_types, - transformed_xstar, - start_time, - settings.time_limit, - work_estimate, - max_work_estimate, - work_limit_reached); - - auto row_combinations = find_mod2_row_combinations( - candidates, max_combination_size, max_row_combinations, &work_estimate, max_work_estimate); - if (work_estimate > max_work_estimate) { work_limit_reached = true; } - scratch_pad_t aggregate_pad(lp.num_cols); - - for (const auto& combination : row_combinations) { - if (toc(start_time) >= settings.time_limit || work_limit_reached) { break; } - - size_t aggregate_input_nz = 0; - for (const i_t candidate_index : combination) { - aggregate_input_nz += candidates[candidate_index].transformed_inequality.size(); - } - const f_t aggregate_work = - (f_t)(4 * aggregate_input_nz + 1) + - (f_t)aggregate_input_nz * std::log2((f_t)aggregate_input_nz + (f_t)1.0); - if (add_work_estimate(aggregate_work, &work_estimate, max_work_estimate, &work_limit_reached)) { - break; - } - - inequality_t aggregate(lp.num_cols); - bool reversible = true; - for (const i_t candidate_index : combination) { - const auto& candidate = candidates[candidate_index]; - aggregate.rhs += candidate.transformed_inequality.rhs; - reversible = reversible && candidate.reversible; - for (i_t k = 0; k < (i_t)candidate.transformed_inequality.size(); ++k) { - aggregate_pad.add_to_pad(candidate.transformed_inequality.index(k), - candidate.transformed_inequality.coeff(k)); - } - } - aggregate_pad.get_pad(aggregate.vector.i, aggregate.vector.x); - aggregate_pad.clear_pad(); - aggregate.sort(); - aggregate.scale((f_t)0.5); - - i_t cuts_added = 0; - mod2_generate_cuts_from_aggregate(complemented_mir, - cut_pool, - lp, - settings, - Arow, - variable_bounds, - var_types, - xstar, - transformed_xstar, - aggregate, - min_violation, - start_time, - work_estimate, - max_work_estimate, - work_limit_reached, - cuts_added); - // if the final inequality is reversable, try the reversed version as well - if (reversible && toc(start_time) < settings.time_limit && !work_limit_reached) { - aggregate.negate(); - mod2_generate_cuts_from_aggregate(complemented_mir, - cut_pool, - lp, - settings, - Arow, - variable_bounds, - var_types, - xstar, - transformed_xstar, - aggregate, - min_violation, - start_time, - work_estimate, - max_work_estimate, - work_limit_reached, - cuts_added); - } - } - return true; -} - -template -bool complemented_mixed_integer_rounding_cut_t::generate_lifted_mixed_binary_cover( - const inequality_t& transformed_inequality, - const std::vector& var_types, - const std::vector& transformed_xstar, - inequality_t& transformed_cut, - f_t& work_estimate, - f_t max_work_estimate) -{ - constexpr f_t tolerance = (f_t)1e-6; - - const f_t estimated_work = - (f_t)(12 * transformed_inequality.size()) + - (f_t)transformed_inequality.size() * std::log2((f_t)transformed_inequality.size() + (f_t)1.0); - if (add_work_estimate(estimated_work, &work_estimate, max_work_estimate)) { return false; } - - inequality_t base = transformed_inequality; - base.negate(); - - std::vector locally_complemented(base.size(), 0); - std::vector solution_value(base.size(), 0.0); - std::vector is_integral(base.size(), 0); - for (i_t k = 0; k < (i_t)base.size(); ++k) { - const i_t j = base.index(k); - f_t aj = base.coeff(k); - if (var_types[j] == variable_type_t::CONTINUOUS) { - solution_value[k] = transformed_xstar[j]; - if (aj > 0.0) { base.vector.x[k] = 0.0; } - continue; - } - - const f_t upper = new_upper(j); - if (upper == inf || std::abs(upper - (f_t)1.0) > tolerance) { return false; } - is_integral[k] = 1; - if (aj < 0.0) { - base.rhs -= aj * upper; - base.vector.x[k] = -aj; - solution_value[k] = upper - transformed_xstar[j]; - locally_complemented[k] = 1; - } else { - solution_value[k] = transformed_xstar[j]; - } - } - - std::vector cover; - cover.reserve(base.size()); - for (i_t k = 0; k < (i_t)base.size(); ++k) { - if (is_integral[k] && base.coeff(k) > tolerance && solution_value[k] > tolerance) { - cover.push_back(k); - } - } - if (cover.empty()) { return false; } - - std::stable_sort( - cover.begin(), cover.end(), lifted_cover_order_t{solution_value, base, tolerance}); - - f_t cover_weight = 0.0; - size_t cover_size = 0; - for (; cover_size < cover.size(); ++cover_size) { - cover_weight += base.coeff(cover[cover_size]); - if (cover_weight - base.rhs > tolerance * std::max((f_t)1.0, std::abs(base.rhs))) { - ++cover_size; - break; - } - } - if (cover_size == 0 || cover_size > cover.size()) { return false; } - cover.resize(cover_size); - - const f_t lambda = cover_weight - base.rhs; - if (lambda <= tolerance) { return false; } - std::sort( - cover.begin(), cover.end(), [&](i_t a, i_t b) { return base.coeff(a) > base.coeff(b); }); - - std::vector prefix(cover.size(), 0.0); - std::vector in_cover(base.size(), 0); - f_t prefix_sum = 0.0; - size_t p = cover.size(); - for (size_t h = 0; h < cover.size(); ++h) { - const i_t k = cover[h]; - in_cover[k] = 1; - if (base.coeff(k) - lambda <= tolerance && p == cover.size()) { p = h; } - if (h < p) { - prefix_sum += base.coeff(k); - prefix[h] = prefix_sum; - } - } - if (p == 0) { return false; } - - transformed_cut = base; - transformed_cut.rhs = -lambda; - for (i_t k = 0; k < (i_t)base.size(); ++k) { - if (!is_integral[k]) { - if (base.coeff(k) >= 0.0) { transformed_cut.vector.x[k] = 0.0; } - continue; - } - if (in_cover[k]) { - transformed_cut.vector.x[k] = std::min(base.coeff(k), lambda); - transformed_cut.rhs += transformed_cut.coeff(k); - } else { - transformed_cut.vector.x[k] = - lifted_cover_coefficient(base.coeff(k), prefix, p, lambda, tolerance); - } - } - - for (i_t k = 0; k < (i_t)transformed_cut.size(); ++k) { - if (!locally_complemented[k]) { continue; } - const i_t j = transformed_cut.index(k); - const f_t coefficient = transformed_cut.coeff(k); - transformed_cut.rhs -= coefficient * new_upper(j); - transformed_cut.vector.x[k] = -coefficient; - } - inequality_t squeezed_cut(transformed_cut.vector.n); - transformed_cut.squeeze(squeezed_cut); - transformed_cut = std::move(squeezed_cut); - transformed_cut.negate(); - return true; -} - -#ifdef DUAL_SIMPLEX_INSTANTIATE_DOUBLE -template bool generate_mod2_zero_half_cuts( - cut_pool_t& cut_pool, - const lp_problem_t& lp, - const simplex_solver_settings_t& settings, - csr_matrix_t& Arow, - const std::vector& new_slacks, - const std::vector& var_types, - const std::vector& xstar, - variable_bounds_t& variable_bounds, - double start_time, - double& work_estimate); - -template bool -complemented_mixed_integer_rounding_cut_t::generate_lifted_mixed_binary_cover( - const inequality_t& transformed_inequality, - const std::vector& var_types, - const std::vector& transformed_xstar, - inequality_t& transformed_cut, - double& work_estimate, - double max_work_estimate); -#endif - -} // namespace cuopt::mathematical_optimization::mip