Closed bknueven closed 4 months ago
Recently when simulating a model, I was ramping a model variable in an attempt to approach the known physical limits (maximum possible performance) of the system. The model became infeasible at a region before the anticipated physical limits in practice.
Knowing this infeasibility was likely caused by something in the model structure, I initialized the model and first optimally solved it at the boundary of the feasible region, then took a small step into the infeasible region and applied the diagnostic tools available.
Based on the approach the model took to attempt to optimally solve the infeasible problem, the tools available (such as display_variables_near_bounds
, display_constraints_with_large_residuals
, and display_extreme_jacobian_entries
) were not pointing to the underlying problem, and instead pointing to very distantly related equations that I could not identify the root cause. Here, compute_infeasibility_explanation
clearly pointed to the problem (in this case an indexed variable was too strictly bounded and a subset of the indices were hitting a bound) that current available methods in the toolbox were not indicating.
I think this would be a valuable tool to include; the main issues I see right now are determining how best to integrate it and finding someone with the time to do so.
Given my understanding of the tool, this is probably something that you would use if and when you get an infeasible solution returned by the solver (i.e. this is not a tool you use when you get an optimal solution but suspect you could get better performance). As such, it probably does not fit directly into the current report_structural_issues
or report_numerical_issues
method, but stands more as a tool to use if you encounter a specific solver failure (in this case infeasible). The question then would be more about how to inform the user that maybe they should try this method, and adding thin wrapper around this tool so it can be called directly from the toolbox.
The one other thing I see is that compute_infeasibility_explanation
uses a logger for all output, which is a bit different from the rest of the diagnostics tools which use a user defined stream. This is probably a minor issue however, although it would mean a little more work.
I think this would be a valuable tool to include; the main issues I see right now are determining how best to integrate it and finding someone with the time to do so.
Given my understanding of the tool, this is probably something that you would use if and when you get an infeasible solution returned by the solver (i.e. this is not a tool you use when you get an optimal solution but suspect you could get better performance). As such, it probably does not fit directly into the current
report_structural_issues
orreport_numerical_issues
method, but stands more as a tool to use if you encounter a specific solver failure (in this case infeasible). The question then would be more about how to inform the user that maybe they should try this method, and adding thin wrapper around this tool so it can be called directly from the toolbox.
At first thought, advanced IDAES users already know to check display_variables_near_bounds
and display_constraints_with_large_residuals
if they suspect their model might be infeasible. Perhaps we could suggest the compute_infeasibility_explanation
as a follow up if neither of these methods readily identify the cause of the potential infeasibility? Admittedly it would be nice to just parse the Ipopt output and look for "converged to a point of local infeasiblity".
The one other thing I see is that
compute_infeasibility_explanation
uses a logger for all output, which is a bit different from the rest of the diagnostics tools which use a user defined stream. This is probably a minor issue however, although it would mean a little more work.
Maybe we should just change this in Pyomo? Given this was just merged into Pyomo/main, I don't foresee an issue with slightly modifying the API.
compute_infeasibility_explanation
as a suggested step if we see constraints with large residuals (although this can blur with poor scaling, but that just means the user needs to fix scaling first).compute_infeasibility_explanation
either way and display the output to the user in whatever form it comes in.I support adding compute_infeasibility_explanation
as a suggestion if we see large residuals.
On the stream front, it seems like we could wrap the pyomo.contrib version as follows (https://stackoverflow.com/questions/1383254/logging-streamhandler-and-standard-streams):
import logging
from pyomo.contrib.iis import mis
def compute_infeasibility_explanation(model, solver, tee=False, tolerance=1e-6, stream=None):
h = logging.StreamHandler(stream)
h.setLevel(logging.INFO)
h.setFormatter(logging.Formatter(fmt=""))
l = logging.getLogger("compute_infeasibility_explanation")
l.setLevel(logging.INFO)
l.addHandler(h)
mis.compute_infeasibility_explanation(model, solver, tee=tee, tolerance=tolerance, logger=l)
One option however is to add compute_infeasibility_explanation as a suggested step if we see constraints with large residuals (although this can blur with poor scaling, but that just means the user needs to fix scaling first).
I support adding
compute_infeasibility_explanation
as a suggestion if we see large residuals.
If we're converged on this I can add the above method to the toolbox, and place the suggestion in the appropriate spot.
Why do you need to set a formatter for the StreamHandler?
Why do you need to set a formatter for the StreamHandler?
Good question -- it isn't needed. Here is the revised implementation, which creates a temporary logger instead.
import logging
from pyomo.contrib.iis import mis
def compute_infeasibility_explanation(model, solver, tee=False, tolerance=1e-6, stream=None):
h = logging.StreamHandler(stream)
h.setLevel(logging.INFO)
l = logging.Logger(name="compute_infeasibility_explanation")
l.setLevel(logging.INFO)
l.addHandler(h)
mis.compute_infeasibility_explanation(model, solver, tee=tee, tolerance=tolerance, logger=l)
I think this plan sounds good - we can finish up any details in the PR.
Sounds good to me too.
The function
compute_infeasibility_explanation
computes sets of mutually infeasible constraints, which can be much more useful for diagnosing infeasibility than the current approach of checking constraint residuals after Ipopt converges to a point of local infeasibility. It was originally developed for WaterTAP (https://github.com/watertap-org/watertap/pull/1289), but a Pyomo developer took the code and put it intopyomo.contrib.iis
(https://github.com/Pyomo/pyomo/pull/3172).It seems like it would be a useful tool to add to the diagnostics toolbox so IDAES users more broadly are aware of is availability. @hunterbarber has successfully leveraged it to identify infeasibility explanations for WaterTAP models.
Below is a basic demo.
output: