0

This is an XY problem, where:

Problem X: My pyomo model, solved using ipopt reports an “optimal” solution. But the values obtained after solving, when substituted into the model, violate 1 or more of the constraints.

Problem Y: How to evaluate Left and Right sides of a pyomo constraint (either inequality or equality), to determine constraint violation? Or how to determine margin of constraint violation? (printing left and right sides of pyomo constraint might be more informative however). Even more useful, would be a method that can take as input an expression, a mapping of variable names, and values (e.g. the reported “optimal” solution by Pyomo), and outputs the expression with the variable evaluated at that value.

The closest I’ve found to this is the “ExpressionReplacementVisitor” Example on the page https://pyomo.readthedocs.io/en/stable/developer_reference/expressions/managing.html#evaluating-expressions. I could envision a method that visits each node in an expression tree, and evaluates it using the given mapping.

Is there an existing method to evaluate an expression using a given mapping? And, is there an easier way to debug Problem “X”?

The model is quite large, or else I’d share some code. Thanks.

makansij
  • 9,303
  • 37
  • 105
  • 183
  • Related: https://stackoverflow.com/questions/51044262/finding-out-reason-of-pyomo-model-infeasibility – makansij Oct 09 '22 at 22:17
  • 1
    Did you see my answer to the question you reference? You can use the technique I showed there to rattle through your constraints and find “negative slack”. … which would be very odd if the solver exited normally. You can always use the generic `display()` command to see the whole mode or individual components with values substituted in like: `model.constraint1.display()` – AirSquid Oct 09 '22 at 22:53
  • Yes I'm seeing that now...trying it now. Thanks @AirSquid. – makansij Oct 09 '22 at 23:18
  • 1
    Awesome. Give a comment back if it doesn’t answer the mail. (And upvote it if it does ;). ) – AirSquid Oct 09 '22 at 23:33
  • I had to change the code in your answer to accommodate indexed constraints, but the general idea works. I'm getting very small negative slack for one of the constraints, so I think it's just a numerical issue. Also, I'm not sure it's necessary to call `set_values` as you did, because the model after solving should already contain those values, but perhaps that was the case for me because the solver returned with termination condition `optimal`. – makansij Oct 10 '22 at 16:28
  • 1
    Agree w/ all. if the delta is decimal dust, it is numerical error of some kind. I only used `set_values` to simulate solving. – AirSquid Oct 10 '22 at 16:49

0 Answers0