Preserve affine GDPopt expressions with mutable-parameter coefficients - #25
Conversation
Use standard representations for GLOA master construction and affine-cut classification so mutable parameter expressions do not hide affine structure. Add regression coverage for objective bounds, enforced constraints, and cut generation. Ref: Pyomo#4036
bernalde
left a comment
bernalde
left a comment
There was a problem hiding this comment.
Reviewed commit 213e420b52eefce58ef3dfae9a8cf1dda6a4f92a against the dedicated submission branch and Pyomo issue Pyomo#4036.
I found no blocking, nonblocking, or question-level code issues. The change consistently uses the standard representation at the three GLOA decision points that previously treated mutable-parameter coefficients as nonlinear, while preserving the existing nonlinear paths. The regression coverage exercises objective preservation, master objective and constraint bounds, and avoidance of MC++ cuts; the full local GDPopt suite passed (79 passed, 35 skipped, 5 deselected). Independent mutation checks made each production classifier fail its corresponding regression.
One current CI failure is external to this diff: the Windows Python 3.10 job crashes with an access violation in the GAMS Python bindings while collecting test_BARON.py. The remaining matrix jobs were still running at review time. This should be rerun before the PR leaves draft.
Finding inventory: 0 Blocking, 0 Nonblocking, 0 Questions.
Merge readiness: the implementation is ready for human review, but the draft/AI-review gate and a clean current-head CI result remain outstanding.
bernalde
left a comment
bernalde
left a comment
There was a problem hiding this comment.
This is a second, independent pass at the same head and supersedes my earlier review of 213e420b. It confirms the three changed sites and adds one scope question the first pass did not raise.
The classifier swap is the right structural fix. For the issue's expression, polynomial_degree() returns None while generate_standard_repn(...).is_linear() is True with linear coefficient 11.810386265582135, so the GLOA master was discarding exact affine structure. Using the standard representation is also the direction this code has already been moving: loa.py:217 has used generate_standard_repn(constr.body).is_linear() for the LOA OA-cut site since 7bac258 ("Switching from polynomial_degree to standard_repn in the OA cut generation so we don't get confused by anything fishy in polynomial_degree"), which retired an explicit # ESJ TODO: This is a risky use of polynomial_degree. This PR converts three more sites of the same kind.
Question: Fixes overstates the scope — four polynomial_degree() routing sites remain
Pyomo#4036 asks for more than the two reproducer cases: "Other GDPopt polynomial_degree() call sites that make solver-routing or linear/nonlinear decisions should be audited for the same failure mode." After this PR, four such sites remain in pyomo/contrib/gdpopt/:
util.py:93,95—solve_continuous_problem()LP-vs-NLP routingsolve_subproblem.py:360—call_appropriate_subproblem_solver()linear-vs-NLP routingdiscrete_problem_initialize.py:292—disjunct_needs_coverin set-covering initializationbranch_and_bound.py:169— deactivation of "nonlinear" constraints in the LBB relaxation
(loa.py:309 asserts on a generated OA cut rather than routing a model, so I read it as out of scope.)
Two of the four reproduce with the issue's own CRF expression.
util.py:93 — a model that is an LP in substance fails outright instead of solving:
m.crf = Expression(expr=m.rate / (1 - (1 + m.rate) ** (-m.years)))
m.c = Constraint(expr=m.x / m.crf <= 9 / m.crf)
m.obj = Objective(expr=m.x / m.crf, sense=maximize)
SolverFactory('gdpopt.gloa').solve(m, mip_solver='highs')At this head that raises RuntimeError: Attempting to use an unavailable solver ... "couenne" — the model is routed to the global NLP solver even though mip_solver='highs' was supplied and every constraint and the objective are affine. Applying the same standard-representation test at that site solves it as an LP and returns 106.29347639023919, the value the issue's table names as correct.
branch_and_bound.py:169 — with the same CRF constraint inside a disjunct, GDPopt-LBB deactivates d_disjuncts[0].constraint[1] from the relaxation and records it in disjunct_to_nonlinear_constraints. That is the same weakened-bound symptom the issue reports, in a different algorithm.
I did not separately exercise solve_subproblem.py:360 or discrete_problem_initialize.py:292.
Either extend the fix to those sites, or change the body's Fixes Pyomo/pyomo#4036 to Refs Pyomo/pyomo#4036 and open a tracked follow-up for the audit before the upstream PR is opened. The branch's commit message already uses Ref: and this PR targets a non-default branch, so nothing auto-closes on merge either way; this is about what the upstream PR body will claim.
Two smaller items are inline.
Verification
At 213e420b, Python 3.13.5, HiGHS 1.14.0:
pyomo/contrib/gdpopt/tests/test_gloa_affine_mutable_param.py: 4 passed. The same four tests run against unmodifiede099f646fail with exactly the issue's numbers —118.10386265582136instead of106.2934763902392, and10.0instead of1.0.- Per-site mutation: restoring each old classifier on its own fails exactly one test and leaves the other three passing —
util.pybreakstest_gloa_does_not_replace_affine_objective,create_oa_subproblems.pybreakstest_gloa_keeps_affine_constraint,gloa.pybreakstest_gloa_does_not_generate_cuts_for_affine_constraint. All three changed sites are independently pinned. pyomo/contrib/gdpopt/tests: 82 passed, 32 skipped, 5 deselected. The skips are the unavailable GAMS/Gurobi/BARON/global-solver lanes.- Parity on the existing corpus: across 2549 constraints and objectives in 20 models under
examples/gdp/, the old and new classifiers disagree zero times, so the change is inert outside the mutable-parameter pattern. - Reverse direction: the only expression class I found where
polynomial_degree()says affine andis_linear()says nonlinear isExpr_ifover a fixed-Varcondition. That path raisesIncompatibleModelError: Highs interface does not support expressions of degree Noneat bothe099f646and this head, so nothing that previously worked is newly dropped. - Cost:
generate_standard_repn(...).is_linear()is faster thanpolynomial_degree()on affine bodies (15 vs 33 µs on a 50-term constraint) and about 4.6x slower on nonlinear ones (5.7 vs 1.2 ms at 1000 terms), which is small next to the MC++ call that immediately follows it in_add_cuts_to_discrete_problem. No performance concern.
CI on this head: 9 checks passing, 4 failing, and all four are unrelated to this change.
lint/style-and-typos— the URL Checker step, on deadminlp.orglinks. Black and typos passed.linux/3.12/conda—test_tee.py::TestTeeStream::test_merged_out_and_err_without_peek, a stream-ordering flake.linux/3.14—test_gdpopt.py::TestGDPopt::test_time_limit, which requires the 8-process model not to converge within 1 s and it converged at 1.03 s. The two classifiers disagree on none of that model's constraints or its objective, so this PR does not move that solve's trajectory.win/3.10/pip— a segfault in the GAMS bindings (cfgmcc.cfgFree) while collectingtest_BARON.py.
The four new tests ran and passed on linux/3.14.
Summary
- Blocking: none.
- Nonblocking: 2, both inline.
- Questions: 1, the
Fixes-versus-Refsscope question above. - Merge-readiness: the code change is correct, minimal, and well covered by tests I was able to break site by site. This stays a COMMENT because the PR is still a draft; the remaining decision is the scope question, then marking it ready.
| Constraint, active=True, descend_into=(Block, Disjunct) | ||
| ): | ||
| if c.body.polynomial_degree() not in (1, 0): | ||
| if not generate_standard_repn(c.body).is_linear(): |
There was a problem hiding this comment.
Nonblocking: nothing at any of the three sites records why polynomial_degree() is not good enough, and the same bare expression now appears in four places (loa.py:217, this line, gloa.py:158, and util.py:116). When this swap was last made, in 7bac258, the # ESJ TODO: This is a risky use of polynomial_degree comment went away with it. A future reader has no signal that "simplifying" this back to the cheaper polynomial_degree() check reintroduces Pyomo#4036.
Suggest one helper in pyomo/contrib/gdpopt/util.py, with the reason in its docstring, that all sites call:
def is_affine(expr):
"""Return True if expr is affine in the model variables.
polynomial_degree() returns None when a coefficient is a nonlinear
function of mutable Params, even though the expression is affine in
every variable; the standard representation evaluates those
coefficients. See Pyomo/pyomo#4036.
"""
return generate_standard_repn(expr).is_linear()util.py imports nothing from gloa.py, loa.py, or create_oa_subproblems.py, and all three already import from it, so there is no cycle. I applied this refactor across the three changed sites locally: test_gloa_affine_mutable_param.py still gives 4 passed and pyomo/contrib/gdpopt/tests still gives 82 passed, 32 skipped, 5 deselected — identical to this head.
It also makes the four remaining routing sites in my review comment a one-line change each if that audit is taken up. A short comment at each site would do instead if you would rather keep the diff minimal.
There was a problem hiding this comment.
Addressed in 8b3f4206f: util.is_affine() now centralizes the union classifier
and documents both the mutable-Param and fixed-condition Expr_if cases. All
GDPopt routing/cut sites, including LOA and the four audit sites, now use it.
| m.choose_x = Disjunction(expr=[[m.x == 1], [m.x == 9]]) | ||
| m.obj = pyo.Objective(expr=m.x / m.crf, sense=pyo.maximize) | ||
|
|
||
| self.assertTrue(generate_standard_repn(m.obj.expr).is_linear()) |
There was a problem hiding this comment.
Nonblocking: this pins the new classifier but not the divergence that motivates the change. assertTrue(generate_standard_repn(...).is_linear()) would still hold if a future Pyomo release taught polynomial_degree() to handle this expression, and at that point the fixture would silently stop exercising the misclassification it exists to guard.
Asserting the other half documents the mechanism and fails loudly if the premise ever moves:
| self.assertTrue(generate_standard_repn(m.obj.expr).is_linear()) | |
| self.assertIsNone(m.obj.expr.polynomial_degree()) | |
| self.assertTrue(generate_standard_repn(m.obj.expr).is_linear()) |
Same for m.limit.body on line 73. I confirmed both hold at this head: m.obj.expr.polynomial_degree() and m.limit.body.polynomial_degree() are each None.
There was a problem hiding this comment.
Addressed in 8b3f4206f: both the objective and constraint regressions now
assert that polynomial_degree() is None before asserting that the standard
representation is linear. The focused module passes all 6 tests.
bernalde
left a comment
bernalde
left a comment
There was a problem hiding this comment.
Correction to my previous review of this head. That review said:
Reverse direction: the only expression class I found where
polynomial_degree()says affine andis_linear()says nonlinear isExpr_ifover a fixed-Varcondition. That path raisesIncompatibleModelError: Highs interface does not support expressions of degree Noneat bothe099f646and this head, so nothing that previously worked is newly dropped.
The second sentence is wrong. I checked that path through HiGHS only, and generalised from one writer. With glpk or cbc — and glpk is the mip_solver the rest of pyomo/contrib/gdpopt/tests uses — the base solves it correctly and this head does not. The expression class is real and the regression is real; details and a verified fix below. The rest of that review stands, including the scope question, which is still open.
Blocking: the change trades one silent bound weakening for another
generate_standard_repn(...).is_linear() and polynomial_degree() are each incomplete, and they fail on disjoint expression classes. Swapping one for the other fixes Pyomo#4036 and opens the mirror-image hole.
Probing 75 expressions built from mutable Params, fixed Vars, external functions and named Expressions, the two classifiers disagree 13 times, and the disagreements fall into exactly two groups:
- 7 where
polynomial_degree()isNoneand the standard representation is linear — the CRF family this PR fixes. - 6 where
polynomial_degree()is1and the standard representation is nonlinear — every one anExpr_ifwhose condition is potentially variable but currently fixed.
The second group is not hypothetical. Expr_ifExpression._compute_polynomial_degree resolves the branch whenever the condition has degree 0 (numeric_expr.py:1455), while _collect_branching_expr only resolves it when the condition is not potentially variable (standard_repn.py:994). A fixed Var is still potentially variable, so the standard representation gives up on an expression that is affine in every variable — the same shape of mistake as Pyomo#4036, in the other direction.
m.x = Var(bounds=(0, 10))
m.sw = Var(bounds=(0, 1), initialize=1)
m.sw.fix(1)
m.limit = Constraint(expr=Expr_if(IF=m.sw >= 0.5, THEN=m.x, ELSE=2 * m.x) <= 1)
m.choose_x = Disjunction(expr=[[m.x >= 0], [m.x == 10]])
m.obj = Objective(expr=m.x, sense=maximize)
SolverFactory('gdpopt.gloa').solve(m, init_algorithm='no_init', iterlim=1, mip_solver=...)Master upper bound after one iteration (exact value is 1.0):
mip_solver |
e099f646 (base) |
213e420b (this head) |
|---|---|---|
glpk |
1.0 |
10.0 |
cbc |
1.0 |
10.0 |
highs |
writer error | 10.0 |
create_oa_subproblems.py:79 now deactivates that constraint from the master. The bound stays valid, but it is exactly the silent weakening Pyomo#4036 is about, and with the default mip_solver it is a regression from working behaviour on this branch. The HiGHS column is a separate pre-existing limitation of that writer and is unchanged by either side.
Fix
Accept an expression that either test can prove affine, and put the reason in one place. Standard-repn first, so the affine path costs what this PR already costs and only genuinely nonlinear bodies pay both:
def is_affine(expr):
"""Return True if expr is affine in the model variables.
Neither available test is complete on its own, and they fail on
disjoint expression classes, so we accept an expression that either
one can prove affine:
- polynomial_degree() returns None when a coefficient is a nonlinear
function of mutable Params, even though the expression is affine in
every variable (Pyomo/pyomo#4036).
- The standard representation does not resolve an Expr_if whose
condition is potentially variable, even when every variable in that
condition is currently fixed and polynomial_degree() resolves it.
"""
return generate_standard_repn(expr).is_linear() or expr.polynomial_degree() in (1, 0)with the three changed sites calling it (util.py:116, create_oa_subproblems.py:79, gloa.py:158). That is the helper from my inline comment on create_oa_subproblems.py, now carrying real logic rather than being a rename.
I applied this and re-ran everything:
- The
Expr_ifmodel above returns1.0again onglpkandcbc, and the pre-existing HiGHS writer error is unchanged. test_gloa_affine_mutable_param.py: 4 passed — the Pyomo#4036 fix is untouched.pyomo/contrib/gdpopt/tests: 82 passed, 32 skipped, 5 deselected, identical to this head.- Over the 75-expression probe the union is exactly the disjunction of the two classifiers, calling 63 affine against 56 and 57. Every one of the 12 it still rejects is genuinely nonlinear (bilinear products,
x**2, or anExpr_ifon a genuinely variable condition). - Cost on a 500-term body: affine
150.4 µsversus149.4 µsfor this PR's current code; nonlinear3381 µsversus2885 µs. In_add_cuts_to_discrete_problemthe extra work is 8–10% of the singlemc()call on the next line, measured with MC++ available.
Regression test
The suite needs a case in this direction, and it cannot use highs, since that writer rejects Expr_if at both base and head. This one passes at e099f646, fails at 213e420b, and passes with the fix:
@unittest.skipUnless(
pyo.SolverFactory('glpk').available(exception_flag=False), 'glpk is not available'
)
class TestGLOAAffineExprIf(unittest.TestCase):
def test_gloa_keeps_affine_expr_if_with_fixed_condition(self):
m = pyo.ConcreteModel()
m.x = pyo.Var(bounds=(0, 10))
m.sw = pyo.Var(bounds=(0, 1), initialize=1)
m.sw.fix(1)
m.limit = pyo.Constraint(
expr=pyo.Expr_if(IF=m.sw >= 0.5, THEN=m.x, ELSE=2 * m.x) <= 1
)
m.choose_x = Disjunction(expr=[[m.x >= 0], [m.x == 10]])
m.obj = pyo.Objective(expr=m.x, sense=pyo.maximize)
# polynomial_degree() resolves the fixed condition; the standard
# representation does not.
self.assertEqual(m.limit.body.polynomial_degree(), 1)
self.assertFalse(generate_standard_repn(m.limit.body).is_linear())
result = pyo.SolverFactory('gdpopt.gloa').solve(
m, init_algorithm='no_init', iterlim=1, mip_solver='glpk'
)
self.assertEqual(
result.solver.termination_condition, TerminationCondition.maxIterations
)
self.assertAlmostEqual(result.problem.upper_bound, 1)That the new file's own solver choice hid this is worth acting on separately: pinning the class to highs means any expression class HiGHS cannot write is invisible to it. The rest of pyomo/contrib/gdpopt/tests uses module-level mip_solver = 'glpk' and nlp_solver = 'ipopt' with skipUnless(SolverFactory(mip_solver).available(), ...), and following that idiom here would have surfaced this.
Base currency
fix/issue-4036-gdpopt-affine-mutable-param is at e099f6463abd3b17b3231abf108b0a46e9640d1a, which is the current tip of Pyomo/pyomo:main — zero commits behind. Nothing to rebase.
Summary
- Blocking: 1, the
Expr_ifregression above. - Nonblocking: 1 still open, the test-assertion comment on
test_gloa_affine_mutable_param.py:38. The helper comment oncreate_oa_subproblems.py:79is superseded by the fix above. - Questions: 1 still open, the
Fixes-versus-Refsscope question from my previous review. - I would not merge this until the blocking issue above is addressed.
| m.component_data_objects(Objective, descend_into=True, active=True) | ||
| ) | ||
| if discrete_obj.polynomial_degree() in (1, 0): | ||
| if generate_standard_repn(discrete_obj.expr).is_linear(): |
There was a problem hiding this comment.
Blocking: this site, and the two others changed here, need the union of both tests rather than a swap. generate_standard_repn(...).is_linear() returns False for an Expr_if whose condition is potentially variable but currently fixed, where polynomial_degree() returns 1 and resolves the branch — so the master silently drops an affine constraint that e099f646 kept. With mip_solver='glpk' the one-iteration master bound on the model in my review comment goes from 1.0 at base to 10.0 here.
is_affine() belongs in util.py next to this call, since all three sites and loa.py:217 want it:
def is_affine(expr):
"""Return True if expr is affine in the model variables.
Neither available test is complete on its own, and they fail on
disjoint expression classes, so we accept an expression that either
one can prove affine:
- polynomial_degree() returns None when a coefficient is a nonlinear
function of mutable Params, even though the expression is affine in
every variable (Pyomo/pyomo#4036).
- The standard representation does not resolve an Expr_if whose
condition is potentially variable, even when every variable in that
condition is currently fixed and polynomial_degree() resolves it.
"""
return generate_standard_repn(expr).is_linear() or expr.polynomial_degree() in (1, 0)then this line becomes if is_affine(discrete_obj.expr):, create_oa_subproblems.py:79 becomes if not is_affine(c.body):, and gloa.py:158 becomes if is_affine(constr.body):. Ordering the standard-repn test first keeps the affine path at the cost this PR already pays (150.4 vs 149.4 µs on a 500-term body); only genuinely nonlinear bodies evaluate both.
I ran this: the Expr_if bound returns to 1.0 on glpk and cbc, test_gloa_affine_mutable_param.py still gives 4 passed, and pyomo/contrib/gdpopt/tests still gives 82 passed, 32 skipped, 5 deselected.
There was a problem hiding this comment.
Addressed in 8b3f4206f: the helper accepts either classifier's affine result,
with standard representation evaluated first. The new GLPK fixed-condition
Expr_if regression reproduced the pre-fix bound of 10 instead of 1 and now
passes with upper bound 1; a standard-representation-only mutation makes it red.
Review feedback addressedPushed commit
This addresses the review-body findings Verification on Python 3.13.5 from the repository
This comment supersedes the PR body's stale four-test, 79-test, four-file, and Correction: the earlier wording overstated the affine-cut test's coverage. A |
Follow-up coverage and review readinessPushed commits
Verification on Python 3.13.5 from the repository
No new inline feedback appeared after the push, and no feedback was |
1 participant
Fixes Pyomo#4036.
This is a staging PR in the contributor fork. It targets a non-default branch,
so merging it will not close the upstream issue. After human review and merge,
fix/issue-4036-gdpopt-affine-mutable-paramis intended to be the source branchfor a separate PR into
Pyomo/pyomo:main.Summary/Motivation:
GDPopt used
polynomial_degree()at several solver-routing andlinear/nonlinear classification sites. Expressions with coefficients that
depend nonlinearly on mutable
Paramobjects can return degreeNoneeventhough they are affine in every model variable. GDPopt could consequently
replace an affine objective, weaken a master or LBB relaxation, send an LP
subproblem to an NLP solver, and perform unnecessary set-cover initialization
or nonlinear cut processing.
The issue's interval objective bounds are valid. The defect is that GDPopt used
those bounds after structurally misclassifying the expression. Neither
polynomial_degree()nor the standard representation is complete alone:fixed-condition
Expr_ifexpressions expose the reverse disagreement. The fixtherefore treats an expression as affine when either classifier can prove it.
Changes proposed in this PR:
util.is_affine()helper that combines the standardrepresentation and polynomial-degree classifiers.
routing, objective replacement, GLOA/LOA cut classification, discrete-master
construction, set-cover initialization, and LBB relaxation construction.
handling, fixed-condition
Expr_if, continuous and post-fixing subproblemroutes, set-cover initialization, cut classification, and LBB preprocessing.
Acceptance criteria from Pyomo#4036:
GLOA master.
for the same classifier divergence.
Expr_ifexpressions remain affine.Verification evidence (Python 3.13.5 from the repository
.venv):python -m pytest pyomo/contrib/gdpopt/tests/test_gloa_affine_mutable_param.py -q -rs:9 passed.
python -m pytest pyomo/contrib/gdpopt/tests -q -rs: 84 passed, 35 skipped,5 deselected. Skips are for unavailable licensed or optional solvers.
commits were rechecked on the current head.
python -m compileall -q pyomo/contrib/gdpopt: passed.Pyomo/pyomo@e099f6463with bounds118.10386265582136instead of106.2934763902392, and10.0instead of1.0.fixed-condition
Expr_iffailures.decisions each failed exactly its intended guard while the other eight
focused tests passed.
Branch Hygiene
fix/issue-4036-gdpopt-affine-mutable-paramate099f6463abd3b17b3231abf108b0a46e9640d1a, currently also the tip ofPyomo/pyomo:main.review/issue-4036-gdpopt-affine-mutable-param, with four commits andeight changed files scoped to this issue.
upstream default branch.
AI-Use Disclosure
or
AI tools contributed to the development of this PR
Review process (select ONE):
Notes for reviewers (optional): This staging PR is ready for human review. Please focus on the union classifier in
util.is_affine()and whether the route-specific tests adequately pin the continuous, subproblem, set-cover, cut-generation, and LBB call sites. The eventual upstream PR should update the review-process selection after the author completes that review.Legal Acknowledgement
By contributing to this software project, I have read the contribution guide and agree to the following terms and conditions for my contribution: