Pyomo / pyomo

An object-oriented algebraic modeling language in Python for structured optimization problems.
https://www.pyomo.org
Other
2.01k stars 517 forks source link

IPOPT solver does not work on PyROS Usage Example #2531

Open makansij opened 2 years ago

makansij commented 2 years ago

Summary

I’ve been experimenting with different solvers to see which ones to use under which robust optimizations. So, I started with the example provided under https://pyomo.readthedocs.io/en/stable/contributed_packages/pyros.html#pyros-usage-example. However, I was unable to get the single-stage version of this problem to work with IPOPT. The code in the example uses the baron solver. However, since there don’t seem to be any integral uncertain parameters, nor integral variables, this problem should be solvable by a non-linear solver, e.g. IPOPT. I tried this example with IPOPT and got the following error:

Steps to reproduce the issue

# === Required import ===
import pyomo.environ as pyo
import pyomo.contrib.pyros as pyros

# === Instantiate the PyROS solver object ===
pyros_solver = pyo.SolverFactory("pyros")

# === Construct the Pyomo model object ===
m = pyo.ConcreteModel()
m.name = "hydro"

# === Define variables ===
m.x1 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x2 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x3 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x4 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x5 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x6 = pyo.Var(within=pyo.Reals,bounds=(150,1500),initialize=150)
m.x7 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x8 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x9 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x10 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x11 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x12 = pyo.Var(within=pyo.Reals,bounds=(0,1000),initialize=0)
m.x13 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x14 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x15 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x16 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x17 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x18 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x19 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x20 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x21 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x22 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x23 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x24 = pyo.Var(within=pyo.Reals,bounds=(0,None),initialize=0)
m.x25 = pyo.Var(within=pyo.Reals,bounds=(100000,100000),initialize=100000)
m.x26 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)
m.x27 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)
m.x28 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)
m.x29 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)
m.x30 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)
m.x31 = pyo.Var(within=pyo.Reals,bounds=(60000,120000),initialize=60000)

# === Define parameters ===
m.set_of_params = pyo.Set(initialize=[0, 1, 2, 3])
nominal_values = {0:82.8*0.0016, 1:4.97, 2:4.97, 3:1800}
m.p = pyo.Param(m.set_of_params, initialize=nominal_values, mutable=True)

# === Specify the objective function ===
m.obj = pyo.Objective(expr=m.p[0]*m.x1**2 + 82.8*8*m.x1 + 82.8*0.0016*m.x2**2 +
                              82.8*82.8*8*m.x2 + 82.8*0.0016*m.x3**2 + 82.8*8*m.x3 +
                              82.8*0.0016*m.x4**2 + 82.8*8*m.x4 + 82.8*0.0016*m.x5**2 +
                              82.8*8*m.x5 + 82.8*0.0016*m.x6**2 + 82.8*8*m.x6 + 248400,
                              sense=pyo.minimize)

# === Specify the constraints ===
m.c2 = pyo.Constraint(expr=-m.x1 - m.x7 + m.x13 + 1200<= 0)
m.c3 = pyo.Constraint(expr=-m.x2 - m.x8 + m.x14 + 1500 <= 0)
m.c4 = pyo.Constraint(expr=-m.x3 - m.x9 + m.x15 + 1100 <= 0)
m.c5 = pyo.Constraint(expr=-m.x4 - m.x10 + m.x16 + m.p[3] <= 0)
m.c6 = pyo.Constraint(expr=-m.x5 - m.x11 + m.x17 + 950 <= 0)
m.c7 = pyo.Constraint(expr=-m.x6 - m.x12 + m.x18 + 1300 <= 0)
m.c8 = pyo.Constraint(expr=12*m.x19 - m.x25 + m.x26 == 24000)
m.c9 = pyo.Constraint(expr=12*m.x20 - m.x26 + m.x27 == 24000)
m.c10 = pyo.Constraint(expr=12*m.x21 - m.x27 + m.x28 == 24000)
m.c11 = pyo.Constraint(expr=12*m.x22 - m.x28 + m.x29 == 24000)
m.c12 = pyo.Constraint(expr=12*m.x23 - m.x29 + m.x30 == 24000)
m.c13 = pyo.Constraint(expr=12*m.x24 - m.x30 + m.x31 == 24000)
m.c14 = pyo.Constraint(expr=-8e-5*m.x7**2 + m.x13 == 0)
m.c15 = pyo.Constraint(expr=-8e-5*m.x8**2 + m.x14 == 0)
m.c16 = pyo.Constraint(expr=-8e-5*m.x9**2 + m.x15 == 0)
m.c17 = pyo.Constraint(expr=-8e-5*m.x10**2 + m.x16 == 0)
m.c18 = pyo.Constraint(expr=-8e-5*m.x11**2 + m.x17 == 0)
m.c19 = pyo.Constraint(expr=-8e-5*m.x12**2 + m.x18 == 0)
m.c20 = pyo.Constraint(expr=-4.97*m.x7 + m.x19 == 330)
m.c21 = pyo.Constraint(expr=-m.p[1]*m.x8 + m.x20 == 330)
m.c22 = pyo.Constraint(expr=-4.97*m.x9 + m.x21 == 330)
m.c23 = pyo.Constraint(expr=-4.97*m.x10 + m.x22 == 330)
m.c24 = pyo.Constraint(expr=-m.p[2]*m.x11 + m.x23 == 330)
m.c25 = pyo.Constraint(expr=-4.97*m.x12 + m.x24 == 330)

# === Specify which parameters are uncertain ===
uncertain_parameters = [m.p] # We can pass IndexedParams this way to PyROS, or as an expanded list per index

# === Define the pertinent data ===
relative_deviation = 0.15
bounds = [(nominal_values[i] - relative_deviation*nominal_values[i],
              nominal_values[i] + relative_deviation*nominal_values[i])
              for i in range(4)]

# === Construct the desirable uncertainty set ===
box_uncertainty_set = pyros.BoxSet(bounds=bounds)

ipopt_solver = pyo.SolverFactory('ipopt')

local_solver = ipopt_solver
global_solver = ipopt_solver

# === Designate which variables correspond to first- and second-stage degrees of freedom ===
first_stage_variables =[m.x1, m.x2, m.x3, m.x4, m.x5, m.x6,
                     m.x19, m.x20, m.x21, m.x22, m.x23, m.x24, m.x31]
second_stage_variables = []
# The remaining variables are implicitly designated to be state variables

# === Call PyROS to solve the robust optimization problem ===
results_1 = pyros_solver.solve(model = m,
                                 first_stage_variables = first_stage_variables,
                                 second_stage_variables = second_stage_variables,
                                 uncertain_params = uncertain_parameters,
                                 uncertainty_set = box_uncertainty_set,
                                 local_solver = local_solver,
                                 global_solver= global_solver,
                                 tee=True,
                                 options = {
                                    "objective_focus": pyros.ObjectiveType.worst_case,
                                    "solve_master_globally": True,
                                    "load_solution":False
                                  })

# === Query results ===
time = results_1.time
iterations = results_1.iterations
termination_condition = results_1.pyros_termination_condition
objective = results_1.final_objective_value
# === Print some results ===
single_stage_final_objective = round(objective,-1)
print("Final objective value: %s" % single_stage_final_objective)
print("PyROS termination condition: %s" % termination_condition)

Error Message

Here's the error message:

===========================================================================================
PyROS: Pyomo Robust Optimization Solver v.1.1.1 
Developed by Natalie M. Isenberg (1), John D. Siirola (2), Chrysanthos E. Gounaris (1) 
(1) Carnegie Mellon University, Department of Chemical Engineering 
(2) Sandia National Laboratories, Center for Computing Research

The developers gratefully acknowledge support from the U.S. Department of Energy's 
Institute for the Design of Advanced Energy Systems (IDAES) 
===========================================================================================
======================================== DISCLAIMER =======================================
PyROS is still under development. 
Please provide feedback and/or report any issues by opening a Pyomo ticket.
===========================================================================================

WARNING: Loading a SolverResults object with a warning status into
    model.name="unknown";
      - termination condition: other
      - message from solver: <undefined>
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Input In [5], in <cell line: 107>()
    103 second_stage_variables = []
    104 # The remaining variables are implicitly designated to be state variables
    105 
    106 # === Call PyROS to solve the robust optimization problem ===
--> 107 results_1 = pyros_solver.solve(model = m,
    108                                  first_stage_variables = first_stage_variables,
    109                                  second_stage_variables = second_stage_variables,
    110                                  uncertain_params = uncertain_parameters,
    111                                  uncertainty_set = box_uncertainty_set,
    112                                  local_solver = local_solver,
    113                                  global_solver= global_solver,
    114                                  tee=True,
    115                                  options = {
    116                                     "objective_focus": pyros.ObjectiveType.worst_case,
    117                                     "solve_master_globally": True,
    118                                     "load_solution":False
    119                                   })
    121 # === Query results ===
    122 time = results_1.time

File ~/opt/anaconda3/lib/python3.9/site-packages/pyomo/contrib/pyros/pyros.py:365, in PyROS.solve(self, model, first_stage_variables, second_stage_variables, uncertain_params, uncertainty_set, local_solver, global_solver, **kwds)
    361 model.add_component(model_data.util_block, util)
    362 # Note:  model.component(model_data.util_block) is util
    363 
    364 # === Validate uncertainty set happens here, requires util block for Cardinality and FactorModel sets
--> 365 validate_uncertainty_set(config=config)
    367 # === Leads to a logger warning here for inactive obj when cloning
    368 model_data.original_model = model

File ~/opt/anaconda3/lib/python3.9/site-packages/pyomo/contrib/pyros/util.py:211, in validate_uncertainty_set(config)
    209 # === Check set validity via boundedness and non-emptiness
    210 if not config.uncertainty_set.is_valid(config=config):
--> 211     raise AttributeError("Invalid uncertainty set detected. Check the uncertainty set object to "
    212                          "ensure non-emptiness and boundedness.")
    214 return

AttributeError: Invalid uncertainty set detected. Check the uncertainty set object to ensure non-emptiness and boundedness.

Information on your system

Pyomo version: 6.4.1 Python version: 3.9 Operating system: MacOS Monterey Solver (if applicable): IPOPT

Additional information

Is this the expected behavior? Or am I failing to see how the problem exceeds the capabilities of a non-linear solver such as IPOPT?

makansij commented 2 years ago

PyROS doesn't properly handle an invalid ipopt.opt file in the working directory, and instead of spitting out a relevant message, it reports AttributeError: Invalid uncertainty set detected....

To reproduce, place an invalid ipopt.op file in your working directory, and run the example above.

shermanjasonaf commented 2 years ago

@makansij Running Python 3.9 with the main branch, I encounter the same exception from my end. I also note that I observe the following error messages through the logger (displayed just before the traceback is printed):

ERROR: Solver (ipopt) returned non-zero return code (155)
ERROR: Solver log: Ipopt 3.13.2: Exception of type: OPTION_INVALID in file
    "IpOptionsList.cpp" at line 726:
     Exception message: Read Option: "example". It is not a valid option.
     Check the list of available options.
    ampl_ipopt.cpp: Error in second Initialize!!!!
ERROR: Solver <pyomo.solvers.plugins.solvers.IPOPT.IPOPT object at
    0x7f3786fd6130> encountered exception attempting to optimize master
    problem in iteration 0

If I run with the Pyomo 6.4.1 release, then the outcome is the same, but the third logging statement does not show, as this logging statement was added after #2515.

Did you encounter a similar series of logger statements, or was only the exception traceback printed to your console?