GDP Hull transform: add handling for cases of constraint functions not well-defined at the origin#3880
GDP Hull transform: add handling for cases of constraint functions not well-defined at the origin#3880sadavis1 wants to merge 46 commits into
Conversation
… on solver particulars. Hopefully also less buggy.
it's probably more likely to be installed than baron
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #3880 +/- ##
========================================
Coverage 90.11% 90.12%
========================================
Files 905 905
Lines 107502 107644 +142
========================================
+ Hits 96878 97016 +138
- Misses 10624 10628 +4
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
bernalde
left a comment
bernalde
left a comment
There was a problem hiding this comment.
Please address the inline comments.
| local_vars[disj].add(var) | ||
| all_local_vars.add(var) | ||
| else: | ||
| vars_to_disaggregate[disjunct].add(var) |
There was a problem hiding this comment.
Blocking: This branch uses the stale outer-loop disjunct variable instead of the current disj. For a variable marked LocalVars while appearing in multiple disjuncts, this only adds the variable to vars_to_disaggregate for the last active disjunct, leaving other disjuncts without a substitute. A simple two-disjunct model with x in both disjuncts and LocalVars entries for both currently fails during transformation with ValueError: No value for uninitialized ScalarVar object x.
Please use the current disjunct and add regression coverage for the generalized-local-vars case. The fix may need to ensure every disjunct in disjuncts_var_appears_in[var] gets a disaggregated variable once the variable is treated as generalized local.
There was a problem hiding this comment.
Good catch. I think all that needs to be done here is to use the correct variable disj instead of disjunct (I really dislike how python keeps that variable in scope after its loop is finished...). I fixed that and added a test for generalized localvars to ensure it behaves as intended.
| f"{e.__class__.__name__} (was prepared for success or a ValueError)." | ||
| ) | ||
| raise | ||
| # Second, try making it well-defined by editing only the regular vars |
There was a problem hiding this comment.
Blocking: _get_well_defined_point() mutates original model variables before calling the heuristic solver, but restores values/fixed states only on the successful path. If the solver is unavailable, errors, or no point is found, the user model is left changed. For example, with a fixed initialized variable, a failed hull transform due to missing gurobipy left x changed from value=5, fixed=True to value=0, fixed=False.
Please wrap the solve/search section in try/finally and restore orig_values and orig_fixed on all exits, including solver exceptions and the GDP_Error path.
There was a problem hiding this comment.
Fixed -- now handled in a finally block, so the original model state is restored even when we exit by throwing
| self.assertEqual(m_pts[m.disjunction][m.a], 0) | ||
| self.assertEqual(m_pts[m.disjunction][m.x], 0) | ||
|
|
||
| @unittest.skipUnless(gurobi_available, "Gurobi is not available") |
There was a problem hiding this comment.
Nonblocking: Most new domain-restriction coverage is skipped when Gurobi is unavailable, including the documented escape hatch for users who pass well_defined_points directly. Please add at least one non-Gurobi test that supplies a ComponentMap manually, for example transforming log(m.x - 1) with {m.disjunction: ComponentMap([(m.x, 2)])}.
There was a problem hiding this comment.
I don't think this is a concern--we have Gurobi in our testing infrastructure.
There was a problem hiding this comment.
It's no trouble so I went ahead and added one anyway
| """, | ||
| ), | ||
| ) | ||
| CONFIG.declare( |
There was a problem hiding this comment.
Nonblocking: The class docstring says the transformation accepts keyword arguments but still lists only perspective_function, EPS, and targets. Since this PR adds public options, please document well_defined_points and well_defined_points_heuristic_solver there or point readers to the generated CONFIG documentation.
There was a problem hiding this comment.
Added those entries to the class doc
…ting needed for correctness here.
|
Merged in the latest changes from main, fixing conflicts, and addressed review comments - thank you David for the review. |
|
(1 character fix; test was throwing the right exception but not for the exact reason I wanted) |
jsiirola
left a comment
jsiirola
left a comment
There was a problem hiding this comment.
A bunch of minor questions. I am comfortable deferring all except for what looks like a bug in the division domain handler.
| orig_values = ComponentMap() | ||
| orig_fixed = ComponentMap() | ||
| for x in itertools.chain(regular_vars, fallback_vars): | ||
| x0_map[x] = 0 # ZeroConstant? |
There was a problem hiding this comment.
I think we have (slowly) moving away from using NumericConstant (including ZeroConstant) anywhere in Pyomo.
| orig_fixed = ComponentMap() | ||
| for x in itertools.chain(regular_vars, fallback_vars): | ||
| x0_map[x] = 0 # ZeroConstant? | ||
| orig_values[x] = value(x, exception=False) |
There was a problem hiding this comment.
If these are all going to be Vars, any reason not to just use x.value?
There was a problem hiding this comment.
Changed to using .value throughout when dealing with Vars
| except ValueError: # ('math domain error') | ||
| pass | ||
| except ZeroDivisionError: | ||
| pass |
There was a problem hiding this comment.
This could be
| except ValueError: # ('math domain error') | |
| pass | |
| except ZeroDivisionError: | |
| pass | |
| except (ValueError, ZeroDivisionError): # ('math domain error') | |
| pass |
| results = self._config.well_defined_points_heuristic_solver.solve( | ||
| test_model, load_solutions=False | ||
| ) | ||
| if results.solver.termination_condition is TerminationCondition.infeasible: | ||
| return False | ||
| if results.solver.status is not SolverStatus.ok: | ||
| raise GDP_Error(f"Unexpected solver status {results.solver.status}.") | ||
| test_model.solutions.load_from(results) |
There was a problem hiding this comment.
I wonder if we should just extend the SolverFactory so that it could be used as a domain validator? It already handles the str -> solver mapping. We would just need to augment it to pass through any solver instances that are instances of registered classes (that way the LegacySolverFactory would admit legacy solvers but not (unwrapped) contrib.solver solvers).
| if results.solver.termination_condition is TerminationCondition.infeasible: | ||
| return False | ||
| if results.solver.status is not SolverStatus.ok: | ||
| raise GDP_Error(f"Unexpected solver status {results.solver.status}.") |
There was a problem hiding this comment.
Would it be more consistent / future-proof to just check_optimal_termination? Or maybe not, because we don't care about optimality -- but then, should we make a peer utility for check_feasible_termination that we can use here?
| # Same LP vs MIP problem as before | ||
| return (arg >= EPS_HEURISTIC,) |
There was a problem hiding this comment.
What if the arg is allowed to be negative? (or by bounds is strictly negative?)
There was a problem hiding this comment.
I think this is a problem if arg is strictly negative. Since it's a heuristic, it's okay so long as it can be positive, but I think this would actually go infeasible if it was non-positive?
There was a problem hiding this comment.
I think the best change here it to just relax the bounds while solving, although it does not eliminate every failure case -- see the discussion below
| if ( | ||
| arg.__class__ in EXPR.native_types | ||
| or not arg.is_potentially_variable() | ||
| or node.name not in _unary_function_handlers |
There was a problem hiding this comment.
I wuold generally prefer "positive inaction" - that is, maintain a list (set) of non-restrictive unary functions, and a map of unary functions that imply domain restrictions, and then raise an exception if the node is not in either set.
There was a problem hiding this comment.
Changed - we now check against a set for all UnaryFunctionExpression. Note that this doesn't cover AbsExpression as it's a subclass, but that's currently the only one
bernalde
left a comment
bernalde
left a comment
There was a problem hiding this comment.
Blocking issues:
- The generated domain constraints for nonzero expressions still force the argument to be positive. This makes valid strictly negative denominators and negative-integer-power bases infeasible during the well-defined-point search.
Nonblocking issues:
- There is one stale comment around the numeric zero used for the origin substitution.
Questions: none.
Tests run:
python -m pytest -q pyomo/gdp/tests/test_hull.py::DomainRestrictionTest pyomo/gdp/tests/test_hull.py::WellDefinedConstraintWalkerTest pyomo/gdp/tests/test_hull.py::TestGeneralizedLocalVars-> 4 passed, 6 skipped.python -m pytest -q pyomo/gdp/tests/test_hull.py-> 151 passed, 10 skipped.python -m black --check --diff pyomo/gdp/plugins/hull.py pyomo/gdp/tests/test_hull.py-> 2 files would be left unchanged.- Direct walker reproduction:
_WellDefinedConstraintGeneratoremits0.0001 <= xfor both1 / xandx**-1whenxhas bounds(-2, -1).
Gurobi is not available in my local environment, so the Gurobi-dependent heuristic tests were skipped locally. Hosted CI is currently green.
I would not merge this until the blocking issues above are addressed.
| return () | ||
| else: | ||
| # return base != 0 | ||
| return (base >= EPS_HEURISTIC,) |
There was a problem hiding this comment.
Blocking: This treats every negative-integer power as if the only acceptable nonzero base is positive. A model with x bounded in [-2, -1] and x**-1 is well-defined, but the helper model gets x >= EPS_HEURISTIC and becomes infeasible before the transformation can find a base point. Please use a shared nonzero-domain helper that chooses base <= -EPS_HEURISTIC when expression bounds rule out the positive side, and add a regression test for a strictly negative base.
There was a problem hiding this comment.
See the discussion below
|
Prepared fixes for the remaining actionable review comments, but I could not push them to the PR head branch. Commits pushed:
Push status:
Main changes prepared:
Tests run:
Comments intentionally not changed in this commit:
Remaining risk:
|
|
This business with positive versus nonzero values is a very good point. Keeping bounds is preferable for a nonlinear solver, and I originally thought it was fine to keep them since it was unlikely to break feasibility of the subproblem (though there are somewhat contrived counterexamples). But since the constraint associated with certain functions was changed from x != 0 to x > 0, it becomes much easier to break feasibility here just by bounding the variable on the wrong side of the number line. There are a few different ways to resolve this, each with some benefits and downsides.
|
6 participants
Fixes # n/a
Summary/Motivation:
When applying the gdp.hull transformation to a model containing a constraint function whose evaluation is not well-defined when every variable is fixed at zero, the perspective functions used by the hull transformation are not well-defined either and lead to errors, especially when the mode is set to FurmanSawayaGrossmann (which is the useful one). This PR slightly alters the mathematical formulation to permit a nonzero base point which can be used in place of the origin, and adds a heuristic to try to magically find one by calling gurobi. When the zero point works, it uses that so Gurobi is never invoked.
Note: this is a minor merge conflict with #3874, if that is merged first I can rebase this
Changes proposed in this PR:
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: