In advance of using PyROS to solve robust optimization problems,
you will need (at least) one local nonlinear programming (NLP) solver
(e.g.,
CONOPT,
IPOPT,
Knitro)
and (at least) one global NLP solver
(e.g.,
BARON,
COUENNE,
SCIP)
installed and licensed on your system.
PyROS can be installed as follows:
Install Pyomo.
PyROS is included in the Pyomo software package, at pyomo/contrib/pyros.
Install NumPy and SciPy with your preferred package manager;
both NumPy and SciPy are required dependencies of PyROS.
You may install NumPy and SciPy with, for example, conda:
conda install numpy scipy
or pip:
(Optional) Test your installation:
install pytest and parameterized
with your preferred package manager (as in the previous step):
pip install pytest parameterized
You may then run the PyROS tests as follows:
python -c 'import os, pytest, pyomo.contrib.pyros as p; pytest.main([os.path.dirname(p.__file__)])'
Some tests involving deterministic NLP solvers may be skipped
if
IPOPT,
BARON,
or
SCIP
is not
pre-installed and licensed on your system.
We now provide a quick overview of how to use PyROS
to solve a robust optimization problem.
Consider the nonconvex deterministic QCQP
\[\begin{array}{clll}
\displaystyle \min_{\substack{x \in [-100, 100], \\ z \in [-100, 100], \\ (y_1, y_2) \in [-100, 100]^2}}
& ~~ x^2 - y_1 z + y_2 & \\
\displaystyle \text{s.t.}
& ~~ xy_1 \geq 150 (q_1 + 1)^2 \\
& ~~ x + y_2^2 \leq 600 \\
& ~~ xz - q_2 y_1 = 2 \\
& ~~ y_1^2 - 2y_2 = 15
\end{array}\]
in which
\(x\) is the sole first-stage variable,
\(z\) is the sole second-stage variable,
\(y_1, y_2\) are the state variables,
and \(q_1, q_2\) are the uncertain parameters.
The uncertain parameters \(q_1, q_2\)
each have a nominal value of 1.
We assume that \(q_1, q_2\)
can independently deviate from their
nominal values by up to \(\pm 10\%\),
so that \((q_1, q_2)\) is constrained in value to the
interval uncertainty set \(\mathcal{Q} = [0.9, 1.1]^2\).
In anticipation of using the PyROS solver and building the deterministic Pyomo
model:
>>> import pyomo.environ as pyo
>>> import pyomo.contrib.pyros as pyros
PyROS can be instantiated through the Pyomo
SolverFactory:
>>> pyros_solver = pyo.SolverFactory("pyros")
We now use PyROS to solve the model to robust optimality
by invoking the solve()
method of the PyROS solver object:
>>> results_1 = pyros_solver.solve(
... # required arguments
... model=m,
... first_stage_variables=first_stage_variables,
... second_stage_variables=second_stage_variables,
... uncertain_params=uncertain_params,
... uncertainty_set=box_uncertainty_set,
... local_solver=local_solver,
... global_solver=global_solver,
... # optional arguments: passed directly to
... # solve to robust optimality
... objective_focus="worst_case",
... solve_master_globally=True,
... )
==============================================================================
PyROS: The Pyomo Robust Optimization Solver...
...
Robust optimal solution identified.
...
All done. Exiting PyROS.
==============================================================================
PyROS, by default, logs to the output console the progress of the optimization
and, upon termination, a summary of the final result.
The summary includes the iteration and solve time requirements,
the final objective function value, and the termination condition.
For further information on the output log,
see the Solver Output Log documentation.
Note
PyROS, like other Pyomo solvers, accepts optional arguments
passed indirectly through the keyword argument options.
This is discussed further in the
Optional Arguments section of the
Solver Interface documentation.
Thus, the PyROS solver invocation in the
preceding code snippet
is equivalent to:
results_1 = pyros_solver.solve(
model=m,
first_stage_variables=first_stage_variables,
second_stage_variables=second_stage_variables,
uncertain_params=uncertain_params,
uncertainty_set=box_uncertainty_set,
local_solver=local_solver,
global_solver=global_solver,
# optional arguments: passed indirectly to
# solve to robust optimality
options={
"objective_focus": "worst_case",
"solve_master_globally": True,
},
)
The PyROS solve() method
returns a results object,
of type ROSolveResults,
that summarizes the outcome of invoking PyROS on a robust optimization problem.
By default, a printout of the results object is included at the end of the solver
output log.
Alternatively, we can display the results object ourselves using:
>>> print(results_1) # output may vary
Termination stats:
Iterations : 3
Solve time (wall s) : 0.917
Final objective value : 9.6616e+03
Termination condition : pyrosTerminationCondition.robust_optimal
We can also query the results object’s individual attributes:
>>> results_1.iterations # total number of iterations
3
>>> results_1.time # total wallclock time; may vary
0.839
>>> results_1.final_objective_value # final objective value; may vary
9661.6...
>>> results_1.pyros_termination_condition # termination condition
<pyrosTerminationCondition.robust_optimal: 1>
Since PyROS has successfully solved our problem,
the final solution has been automatically loaded to the model.
We can inspect the resulting state of the model
by invoking, for example, m.display() or m.pprint().
For a general discussion of the PyROS solver outputs,
see the
Overview of Outputs section of the
Solver Interface documentation.
PyROS uses polynomial decision rules to approximate the adjustability
of the second-stage variables to the uncertain parameters.
The degree of the decision rule polynomials is
specified through the optional keyword argument
decision_rule_order to the PyROS
solve() method.
By default, decision_rule_order is set to 0,
so that static decision rules are used.
Increasing the decision rule order
may yield a solution with better quality:
>>> results_2 = pyros_solver.solve(
... model=m,
... first_stage_variables=first_stage_variables,
... second_stage_variables=second_stage_variables,
... uncertain_params=uncertain_params,
... uncertainty_set=box_uncertainty_set,
... local_solver=local_solver,
... global_solver=global_solver,
... objective_focus="worst_case",
... solve_master_globally=True,
... decision_rule_order=1, # use affine decision rules
... )
==============================================================================
PyROS: The Pyomo Robust Optimization Solver...
...
Robust optimal solution identified.
...
All done. Exiting PyROS.
==============================================================================
Inspecting the results:
>>> print(results_2) # output may vary
Termination stats:
Iterations : 5
Solve time (wall s) : 1.956
Final objective value : 6.5403e+03
Termination condition : pyrosTerminationCondition.robust_optimal
Notice that when we switch from optimizing over static decision rules
to optimizing over affine decision rules,
there is a ~32% decrease in the final objective
value, albeit at some additional computational expense.
In conjunction with standard Pyomo control flow tools,
PyROS facilitates an analysis of the “price of robustness”,
which we define to be the increase in the robust optimal objective value
relative to the deterministically optimal objective value.
Let us, for example, consider optimizing robustly against
an interval uncertainty set \([1 - p, 1 + p]^2\),
where \(p\) is the half-length of the interval.
We can optimize against intervals of increasing half-length \(p\)
by iterating over select values for \(p\) in a for loop,
and in each iteration, solving a robust optimization problem
subject to a corresponding
BoxSet instance:
>>> results_dict = dict()
>>> for half_length in [0.0, 0.1, 0.2, 0.3, 0.4]:
... print(f"Solving problem for {half_length=}:")
... results_dict[half_length] = pyros_solver.solve(
... model=m,
... first_stage_variables=first_stage_variables,
... second_stage_variables=second_stage_variables,
... uncertain_params=uncertain_params,
... uncertainty_set=pyros.BoxSet(
... bounds=[(1 - half_length, 1 + half_length)] * 2
... ),
... local_solver=local_solver,
... global_solver=global_solver,
... objective_focus="worst_case",
... solve_master_globally=True,
... decision_rule_order=1,
... )
...
Solving problem for half_length=0.0:
...
Solving problem for half_length=0.1:
...
Solving problem for half_length=0.2:
...
Solving problem for half_length=0.3:
...
Solving problem for half_length=0.4:
...
All done. Exiting PyROS.
==============================================================================
Using the dict populated in the loop,
and the
previously evaluated deterministically optimal
objective value,
we can print a tabular summary of the results:
>>> for idx, (half_length, res) in enumerate(results_dict.items()):
... if idx == 0:
... # print table header
... print("=" * 71)
... print(
... f"{'Half-Length':15s}"
... f"{'Termination Cond.':21s}"
... f"{'Objective Value':18s}"
... f"{'Price of Rob. (%)':17s}"
... )
... print("-" * 71)
... # print table row
... obj_value, percent_obj_increase = float("nan"), float("nan")
... is_robust_optimal = (
... res.pyros_termination_condition
... == pyros.pyrosTerminationCondition.robust_optimal
... )
... is_robust_infeasible = (
... res.pyros_termination_condition
... == pyros.pyrosTerminationCondition.robust_infeasible
... )
... if is_robust_optimal:
... # compute the price of robustness
... obj_value = res.final_objective_value
... price_of_robustness = (
... (res.final_objective_value - deterministic_obj)
... / deterministic_obj
... )
... elif is_robust_infeasible:
... # infinite objective
... obj_value, price_of_robustness = float("inf"), float("inf")
... print(
... f"{half_length:<15.1f}"
... f"{res.pyros_termination_condition.name:21s}"
... f"{obj_value:<18.2f}"
... f"{100 * price_of_robustness:<.2f}"
... )
... print("-" * 71)
...
=======================================================================
Half-Length Termination Cond. Objective Value Price of Rob. (%)
-----------------------------------------------------------------------
0.0 robust_optimal 5407.94 -0.00
-----------------------------------------------------------------------
0.1 robust_optimal 6540.31 20.94
-----------------------------------------------------------------------
0.2 robust_optimal 7838.50 44.94
-----------------------------------------------------------------------
0.3 robust_optimal 9316.88 72.28
-----------------------------------------------------------------------
0.4 robust_infeasible inf inf
-----------------------------------------------------------------------
The table shows the response of the PyROS termination condition,
final objective value, and price of robustness
to the half-length \(p\).
Observe that:
The optimal objective value for the interval of half-length
\(p=0\) is equal to the optimal deterministic objective value
The objective value (and thus, the price of robustness)
increases with the half-length
For large enough half-length (\(p=0.4\)), the problem
is robust infeasible
Therefore, this example clearly illustrates the potential
impact of the uncertainty set size on the robust optimal
objective function value
and the ease of analyzing the price of robustness
for a given optimization problem under uncertainty.