From b8b1956d27f2c2f7e41ab3bedaf968c009603069 Mon Sep 17 00:00:00 2001 From: Dave Woodruff Date: Mon, 17 Aug 2026 17:14:35 -0700 Subject: [PATCH 1/2] Clear outer_bound when a solve produces none, instead of keeping a stale one A solve that reports no outer bound left the subproblem's previous bound in place. The comment defending this said the stale value "is still a valid outer bound (any Lagrangian dual value bounds the original problem)". That is true per scenario and irrelevant to the sum Ebound actually forms. Ebound computes Sum_s p_s * outer_bound_s, and a Lagrangian bound is valid as a *sum* only when every scenario uses weights from one generation satisfying Sum_s p_s W_s = 0: Sum_s p_s L_s(W'_s) <= Sum_s p_s [f_s(x*) + W'_s^T x*] = OPT + (Sum_s p_s W'_s)^T xbar* Mixing one scenario's stale bound with the others' fresh ones leaves the trailing term in place, so the result is not an outer bound at all. Nothing downstream can catch it. Ebound's collective missing-bound check tests `outer_bound is None`, and the fallback deliberately leaves a number behind, so the check passes. It catches "never computed"; the fallback converts "not computed this time" into "computed earlier", which is exactly the state it cannot distinguish. The hub then latches the best outer bound it has seen (spcommunicator.py, `new > old` when minimizing), so one too-good value is permanent and no later honest bound corrects it. Reachability is narrow but real. no_outer_bound_results screens only infeasible, infeasibleOrUnbounded, unbounded, and a missing results object, so a subproblem has to go bad at iteration k having been fine earlier. Two routes: lagrangian() calls receive_nonant_bounds() before every solve, and tightened nonant bounds from another spoke change the feasible region, so a subproblem genuinely can newly become infeasible; and an unbounded subproblem has Lagrangian value -inf, so a retained finite number over-claims regardless of any weight mixing. The same exposure existed on the ordinary path: the not_good_enough_results branch also left outer_bound untouched, just without a comment claiming that was fine. Both branches now set it to None, which is what _set_initial_bounds already documents as the "not computed" state; Ebound propagates the None and the bound spokes decline to send, as they already do for a first solve that produces no bound. Callers were checked: the Lagrangian spoke, the subgradient hub, and PHBase's trivial_bound all already guard against None. PHBase.post_solve_bound returns the value directly and can now return None, so its docstring says so. Tests are added to test_outer_bound_only.py (already wired into run_coverage.bash and CI), stubbing the failing solve so no solver is needed. All four fail on the unfixed code and pass with the change. Noticed while writing the design for #838, but independent of it: ipopt reports Lower_bound = -inf rather than None, so it takes the other branch and never reaches this path. Co-Authored-By: Claude Opus 5 --- mpisppy/phbase.py | 6 +- mpisppy/spopt.py | 22 ++++-- mpisppy/tests/test_outer_bound_only.py | 97 ++++++++++++++++++++++++++ 3 files changed, 118 insertions(+), 7 deletions(-) diff --git a/mpisppy/phbase.py b/mpisppy/phbase.py index 59a942120..e67696b28 100644 --- a/mpisppy/phbase.py +++ b/mpisppy/phbase.py @@ -547,8 +547,10 @@ def post_solve_bound(self, solver_options=None, verbose=False): If True, displays verbose output. Default False. Returns: - float: - An outer bound on the optimal objective function value. + float or None: + An outer bound on the optimal objective function value, or + None if any subproblem solve produced no bound (the + expectation cannot be formed if a scenario is missing one). Note: This function overwrites current variable values. This is only diff --git a/mpisppy/spopt.py b/mpisppy/spopt.py index f81eb475a..60a0efbdf 100644 --- a/mpisppy/spopt.py +++ b/mpisppy/spopt.py @@ -381,11 +381,19 @@ def _vb(msg): raise if outer_bound is None: - # Leave outer_bound at its previous value rather than - # publish whatever the solver left in the results object. - # The stale bound is still a valid outer bound (any - # Lagrangian dual value bounds the original problem), - # which is what the pre-outer_bound_only path did here. + # This solve produced no bound, so say so. Keeping the + # previous value would be wrong, even though each stale + # value was a valid bound when it was computed: Ebound + # sums p_s * outer_bound_s across scenarios, and a + # Lagrangian bound is only valid as a *sum* when every + # scenario uses weights from one generation satisfying + # sum_s p_s W_s = 0. Mixing a stale bound with fresh ones + # gives sum_s p_s L_s(W'_s) <= OPT + (sum_s p_s W'_s)^T xbar*, + # whose trailing term does not vanish -- so the result is + # not an outer bound at all, and nothing downstream can + # tell. None instead lets Ebound's missing-bound check fire + # and the spoke decline to send. + s._mpisppy_data.outer_bound = None if gripe: print (f"[{self._get_cylinder_name()}] No outer bound for scenario {s.name}") if results is not None: @@ -399,6 +407,10 @@ def _vb(msg): elif sputils.not_good_enough_results(results): s._mpisppy_data.solution_available = False + # Same reasoning as the outer_bound_only branch above: a failed + # solve computed no bound, and leaving the previous one in place + # lets Ebound form a mixed-generation sum that is not a bound. + s._mpisppy_data.outer_bound = None if gripe: print (f"[{self._get_cylinder_name()}] Solve failed for scenario {s.name}") diff --git a/mpisppy/tests/test_outer_bound_only.py b/mpisppy/tests/test_outer_bound_only.py index 41cec5d2e..59bf1cd96 100644 --- a/mpisppy/tests/test_outer_bound_only.py +++ b/mpisppy/tests/test_outer_bound_only.py @@ -19,6 +19,7 @@ import unittest import numpy as np +from pyomo.opt import SolverResults, SolverStatus, TerminationCondition import mpisppy.opt.ph import mpisppy.phbase @@ -190,5 +191,101 @@ def test_normal_solve_still_loads_solution(self): self.assertTrue(np.isfinite(s._mpisppy_data.inner_bound)) + +class _NoBoundPlugin: + """A solver plugin that comes back with no usable bound. + + Only what solve_one touches: an options mapping and a solve() returning a + results object whose termination condition is one of the outcomes + no_outer_bound_results screens out. Not a persistent solver, so solve_one + takes the ordinary path. + """ + + def __init__(self, termination_condition): + self.options = {} + self._tc = termination_condition + + def solve(self, s, **kwargs): + results = SolverResults() + results.solver.status = SolverStatus.warning + results.solver.termination_condition = self._tc + return results + + +class TestStaleBoundNotRetained(unittest.TestCase): + """A solve that produces no bound must clear the subproblem's bound. + + Keeping the previous value looks harmless -- each stale value really was a + valid outer bound for the weights it was computed with -- but Ebound sums + p_s * outer_bound_s across scenarios, and a Lagrangian bound is valid as a + *sum* only when every scenario uses weights from a single generation + satisfying sum_s p_s W_s = 0: + + sum_s p_s L_s(W'_s) <= OPT + (sum_s p_s W'_s)^T xbar* + + Mixing one scenario's stale bound with the others' fresh ones leaves that + trailing term in place, so the result is not an outer bound at all. Nothing + downstream can detect it: Ebound's missing-bound check looks for None, and a + retained number is not None. The hub then latches the best outer bound it + has seen, so a single too-good value is never corrected by later honest + ones. + + No solver needed: the failing solve is stubbed. + """ + + def _ph_with_prior_bounds(self, value=-100.0): + """A PH object whose subproblems already carry bounds from an earlier + iteration's weights.""" + ph = _make_ph() + for s in ph.local_scenarios.values(): + s._mpisppy_data.outer_bound = value + return ph + + def test_bound_only_infeasible_clears_the_bound(self): + ph = self._ph_with_prior_bounds() + name = list(ph.local_scenarios)[0] + s = ph.local_scenarios[name] + s._solver_plugin = _NoBoundPlugin(TerminationCondition.infeasible) + ph.solve_one(None, name, s, gripe=False, need_solution=False, + outer_bound_only=True) + self.assertIsNone(s._mpisppy_data.outer_bound) + + def test_bound_only_unbounded_clears_the_bound(self): + # An unbounded subproblem has Lagrangian value -inf, so a retained + # finite number over-claims regardless of any weight mixing. + ph = self._ph_with_prior_bounds() + name = list(ph.local_scenarios)[0] + s = ph.local_scenarios[name] + s._solver_plugin = _NoBoundPlugin(TerminationCondition.unbounded) + ph.solve_one(None, name, s, gripe=False, need_solution=False, + outer_bound_only=True) + self.assertIsNone(s._mpisppy_data.outer_bound) + + def test_ebound_declines_rather_than_mixing_generations(self): + # The end-to-end shape of the bug: one scenario keeps a bound from the + # previous weights while the other two get fresh ones. Ebound must + # refuse, not return a finite number that is not a bound. + ph = self._ph_with_prior_bounds() + names = list(ph.local_scenarios) + stale = ph.local_scenarios[names[0]] + stale._solver_plugin = _NoBoundPlugin(TerminationCondition.infeasible) + ph.solve_one(None, names[0], stale, gripe=False, need_solution=False, + outer_bound_only=True) + for n in names[1:]: + ph.local_scenarios[n]._mpisppy_data.outer_bound = -50.0 + self.assertIsNone(ph.Ebound()) + + def test_failed_ordinary_solve_also_clears_the_bound(self): + # The same exposure exists off the outer_bound_only path: the + # not_good_enough_results branch used to leave outer_bound untouched. + ph = self._ph_with_prior_bounds() + name = list(ph.local_scenarios)[0] + s = ph.local_scenarios[name] + s._solver_plugin = _NoBoundPlugin(TerminationCondition.infeasible) + ph.solve_one(None, name, s, gripe=False, need_solution=False) + self.assertIsNone(s._mpisppy_data.outer_bound) + self.assertFalse(s._mpisppy_data.solution_available) + + if __name__ == "__main__": unittest.main() From ffffe12022380132e277a62aa84e81542d5d429d Mon Sep 17 00:00:00 2001 From: Dave Woodruff Date: Mon, 17 Aug 2026 19:20:09 -0700 Subject: [PATCH 2/2] Column generation: report no bound rather than summing a missing reduced cost Clearing outer_bound on a failed solve made a latent path in column generation reachable, and CI found it: File "mpisppy/cgbase.py", line 567, in add_columns_to_mp_from_results sum_redcosts+=red_cost TypeError: unsupported operand type(s) for +=: 'int' and 'NoneType' build_columns_from_subproblem_solutions already leaves red_cost at None when a subproblem produced no outer bound -- it tests the value for finiteness for exactly that reason, and its comment says so -- but the consumer added it unconditionally. Nothing had reached that path before because outer_bound was never None there; now it can be. Skipping the missing term would be worse than the crash. The caller computes rmp_obj_val + sum(reduced costs) and treats it as an outer bound, so dropping a term overstates the bound -- the same unsoundness this branch is about, one subsystem over. So add_columns_to_mp_from_results returns None when any reduced cost is missing, and both consumers (cgbase.CG and opt/dualcg) leave the incumbent best bound and the convergence metric untouched for that iteration rather than computing from a partial sum. The columns are still added; they are useful whether or not a bound came with them. While guarding dualcg, found that update_dual_center compared outer_bound_candidate against the stored center before testing it for None. The old ordering only avoided a TypeError while the center was still unset, so a None candidate would have raised there once a center existed. The None check now comes first. Tested deterministically rather than through the solver: the CI trigger is a size-limited cplex build whose subproblem solves fail where a full license succeeds, which does not reproduce locally. add_columns_to_mp_from_results only calls self.add_column_for_scenario, so a recorder stands in for the CG object and the None path is exercised directly -- no solver, no MPI. Both new assertions fail against the unfixed code with the TypeError above. Co-Authored-By: Claude Opus 5 --- mpisppy/cgbase.py | 56 +++++++++++++++++++++++--------- mpisppy/opt/dualcg.py | 39 ++++++++++++++-------- mpisppy/tests/test_cg_main.py | 61 +++++++++++++++++++++++++++++++++++ 3 files changed, 128 insertions(+), 28 deletions(-) diff --git a/mpisppy/cgbase.py b/mpisppy/cgbase.py index 613edc123..768968fdd 100644 --- a/mpisppy/cgbase.py +++ b/mpisppy/cgbase.py @@ -555,21 +555,37 @@ def add_columns_to_mp_from_results(self, all_results): reduced cost, scenario cost, and variable values. Returns: - float: The sum of reduced costs + float or None: The sum of reduced costs, or None if any subproblem + reported no reduced cost. + + build_columns_from_subproblem_solutions leaves red_cost at None when a + subproblem's solve produced no outer bound. Skipping such a term and + returning the sum of the rest would overstate the bound -- the caller + adds it to rmp_obj_val and treats the result as an outer bound -- so the + whole sum is reported unavailable instead. The columns themselves are + still added: they are useful regardless of whether a bound came with + them. """ sum_redcosts=0 + bound_available=True for rank_results in all_results: if rank_results is not None: if isinstance(rank_results, list): for result in rank_results: sname, red_cost, scen_cost, xvec = result self.add_column_for_scenario(sname, scen_cost, xvec) - sum_redcosts+=red_cost + if red_cost is None: + bound_available=False + else: + sum_redcosts+=red_cost else: sname, red_cost, scen_cost, xvec = rank_results self.add_column_for_scenario(sname, scen_cost, xvec) - sum_redcosts+=red_cost - return sum_redcosts + if red_cost is None: + bound_available=False + else: + sum_redcosts+=red_cost + return sum_redcosts if bound_available else None def _build_columns_from_xhat_list(self, xhat_list): """ @@ -674,22 +690,32 @@ def iterk_loop(self): sum_redcosts=self.add_columns_to_mp_from_results(all_results) self.add_columns_to_mp_from_results(all_results_xhat_recent) self.add_columns_to_mp_from_results(all_results_xfeas) - self.outer_bound_candidate = self.rmp_obj_val + sum_redcosts + # A None sum means some subproblem produced no reduced cost, so + # rmp_obj_val + sum(reduced costs) is not an outer bound this + # iteration. Leave the incumbent best bound and the convergence + # metric where they are rather than computing from a partial + # sum; the next iteration may well have every term. + if sum_redcosts is None: + self.outer_bound_candidate = None + else: + self.outer_bound_candidate = self.rmp_obj_val + sum_redcosts # Update bounds and convergence metric if self.problem_is_lp: self.best_solution_obj_val= self.rmp_obj_val - if self.best_bound_obj_val is None: - self.best_bound_obj_val = self.outer_bound_candidate - elif self.is_minimizing: - self.best_bound_obj_val = max(self.outer_bound_candidate, self.best_bound_obj_val) - else: - self.best_bound_obj_val = min(self.outer_bound_candidate, self.best_bound_obj_val) - if self.is_minimizing: - self.conv=(self.rmp_obj_val-self.best_bound_obj_val)/abs(self.best_bound_obj_val) - else: - self.conv=(self.best_bound_obj_val-self.rmp_obj_val)/abs(self.best_bound_obj_val) + if self.outer_bound_candidate is not None: + if self.best_bound_obj_val is None: + self.best_bound_obj_val = self.outer_bound_candidate + elif self.is_minimizing: + self.best_bound_obj_val = max(self.outer_bound_candidate, self.best_bound_obj_val) + else: + self.best_bound_obj_val = min(self.outer_bound_candidate, self.best_bound_obj_val) + if self.best_bound_obj_val is not None: + if self.is_minimizing: + self.conv=(self.rmp_obj_val-self.best_bound_obj_val)/abs(self.best_bound_obj_val) + else: + self.conv=(self.best_bound_obj_val-self.rmp_obj_val)/abs(self.best_bound_obj_val) if dprogress: print("") diff --git a/mpisppy/opt/dualcg.py b/mpisppy/opt/dualcg.py index de403b786..a29370f3a 100644 --- a/mpisppy/opt/dualcg.py +++ b/mpisppy/opt/dualcg.py @@ -130,6 +130,13 @@ def update_dual_center(self): Update the regularization center if the outer bound candidate improved. """ + # No bound this iteration means there is nothing to compare against and + # nothing to recenter on. Checked first: the comparison below would + # otherwise raise on None once a center has been established, which the + # old ordering only avoided while the center was still unset. + if self.outer_bound_candidate is None: + return + if self.dual_center_outer_bound is not None: if self.is_minimizing: improved = self.outer_bound_candidate >= self.dual_center_outer_bound @@ -137,8 +144,6 @@ def update_dual_center(self): improved = self.outer_bound_candidate <= self.dual_center_outer_bound if not improved: return - elif self.outer_bound_candidate is None: - return m = self.mp @@ -241,23 +246,31 @@ def iterk_loop(self): sum_redcosts=self.add_columns_to_mp_from_results(all_results) self.add_columns_to_mp_from_results(all_results_xhat_recent) self.add_columns_to_mp_from_results(all_results_xfeas) - self.outer_bound_candidate = self.rmp_obj_val + sum_redcosts + # See cgbase.add_columns_to_mp_from_results: None means some + # subproblem reported no reduced cost, so there is no outer + # bound to form from a partial sum this iteration. + if sum_redcosts is None: + self.outer_bound_candidate = None + else: + self.outer_bound_candidate = self.rmp_obj_val + sum_redcosts if not hasattr(self, "dual_center_outer_bound"): self.dual_center_outer_bound = None self.update_dual_center() # Update outer bound and convergence metric. - if self.best_bound_obj_val is None: - self.best_bound_obj_val = self.outer_bound_candidate - elif self.is_minimizing: - self.best_bound_obj_val = max(self.outer_bound_candidate, self.best_bound_obj_val) - else: - self.best_bound_obj_val = min(self.outer_bound_candidate, self.best_bound_obj_val) - if self.is_minimizing: - self.conv=(self.rmp_obj_val-self.best_bound_obj_val)/abs(self.best_bound_obj_val) - else: - self.conv=(self.best_bound_obj_val-self.rmp_obj_val)/abs(self.best_bound_obj_val) + if self.outer_bound_candidate is not None: + if self.best_bound_obj_val is None: + self.best_bound_obj_val = self.outer_bound_candidate + elif self.is_minimizing: + self.best_bound_obj_val = max(self.outer_bound_candidate, self.best_bound_obj_val) + else: + self.best_bound_obj_val = min(self.outer_bound_candidate, self.best_bound_obj_val) + if self.best_bound_obj_val is not None: + if self.is_minimizing: + self.conv=(self.rmp_obj_val-self.best_bound_obj_val)/abs(self.best_bound_obj_val) + else: + self.conv=(self.best_bound_obj_val-self.rmp_obj_val)/abs(self.best_bound_obj_val) if dprogress: print("") diff --git a/mpisppy/tests/test_cg_main.py b/mpisppy/tests/test_cg_main.py index cdb499fd6..12dfeb845 100644 --- a/mpisppy/tests/test_cg_main.py +++ b/mpisppy/tests/test_cg_main.py @@ -19,6 +19,7 @@ from mpisppy.tests.examples.sizes.sizes import scenario_creator as sizes_creator, \ scenario_denouement as sizes_denouement from mpisppy.tests.utils import get_solver, limit_solver_threads +from mpisppy.cgbase import CGBase import mpisppy.MPI as mpi solver_available, solver_name, persistent_available, persistent_solver_name = get_solver(persistent_OK=False) @@ -244,5 +245,65 @@ def test_sizes_iter0_creates_initial_columns(self): self.assertGreater(cg.next_col[sname], 0) + +class TestRedCostSumming(unittest.TestCase): + """add_columns_to_mp_from_results with a subproblem that reported no bound. + + build_columns_from_subproblem_solutions already leaves red_cost at None + when a subproblem's solve produced no outer bound -- it tests the value for + finiteness precisely because that can happen. The consumer used to add it + unconditionally, which raised + + TypeError: unsupported operand type(s) for +=: 'int' and 'NoneType' + + the moment the case actually occurred (seen in CI, where the size-limited + cplex build fails solves that succeed on a full license). + + Silently skipping the missing term would be worse than the crash: the caller + adds the sum to rmp_obj_val and treats the result as an outer bound, so + dropping a term overstates it. The sum reports itself unavailable instead. + + No solver needed -- the method only calls self.add_column_for_scenario, so a + recorder stands in for the CG object. + """ + + class _Recorder: + def __init__(self): + self.added = [] + + def add_column_for_scenario(self, sname, scen_cost, xvec): + self.added.append(sname) + + def _sum(self, all_results): + rec = self._Recorder() + total = CGBase.add_columns_to_mp_from_results(rec, all_results) + return total, rec.added + + def test_all_bounds_present_sums_them(self): + results = [[("s0", 1.5, 10.0, {}), ("s1", 2.5, 20.0, {})]] + total, added = self._sum(results) + self.assertEqual(total, 4.0) + self.assertEqual(added, ["s0", "s1"]) + + def test_a_missing_bound_makes_the_sum_unavailable(self): + results = [[("s0", 1.5, 10.0, {}), ("s1", None, 20.0, {})]] + total, added = self._sum(results) + self.assertIsNone(total) + # The columns are still added: they are useful whether or not a bound + # came with them. + self.assertEqual(added, ["s0", "s1"]) + + def test_single_tuple_per_rank_is_handled_too(self): + # Ranks may report one tuple rather than a list of them. + self.assertEqual(self._sum([("s0", 3.0, 10.0, {})])[0], 3.0) + self.assertIsNone(self._sum([("s0", None, 10.0, {})])[0]) + + def test_none_rank_results_are_skipped(self): + total, added = self._sum([None, [("s0", 1.0, 10.0, {})]]) + self.assertEqual(total, 1.0) + self.assertEqual(added, ["s0"]) + + + if __name__ == '__main__': unittest.main()