PyROS Solver
PyROS (Pyomo Robust Optimization Solver) is a Pyomobased metasolver for nonconvex, twostage adjustable robust optimization problems.
It was developed by Natalie M. Isenberg, Jason A. F. Sherman, and Chrysanthos E. Gounaris of Carnegie Mellon University, in collaboration with John D. Siirola of Sandia National Labs. The developers gratefully acknowledge support from the U.S. Department of Energy’s Institute for the Design of Advanced Energy Systems (IDAES).
Methodology Overview
Below is an overview of the type of optimization models PyROS can accommodate.
PyROS is suitable for optimization models of continuous variables that may feature nonlinearities (including nonconvexities) in both the variables and uncertain parameters.
PyROS can handle equality constraints defining state variables, including implicit state variables that cannot be eliminated via reformulation.
PyROS allows for twostage optimization problems that may feature both firststage and secondstage degrees of freedom.
PyROS is designed to operate on deterministic models of the general form
where:
\(x \in \mathcal{X}\) are the “design” variables (i.e., firststage degrees of freedom), where \(\mathcal{X} \subseteq \mathbb{R}^{n_x}\) is the feasible space defined by the model constraints (including variable bounds specifications) referencing \(x\) only.
\(z \in \mathbb{R}^{n_z}\) are the “control” variables (i.e., secondstage degrees of freedom)
\(y \in \mathbb{R}^{n_y}\) are the “state” variables
\(q \in \mathbb{R}^{n_q}\) is the vector of model parameters considered uncertain, and \(q^{\text{nom}}\) is the vector of nominal values associated with those.
\(f_1\left(x\right)\) are the terms of the objective function that depend only on design variables
\(f_2\left(x, z, y; q\right)\) are the terms of the objective function that depend on all variables and the uncertain parameters
\(g_i\left(x, z, y; q\right)\) is the \(i^\text{th}\) inequality constraint function in set \(\mathcal{I}\) (see Note)
\(h_j\left(x, z, y; q\right)\) is the \(j^\text{th}\) equality constraint function in set \(\mathcal{J}\) (see Note)
Note
PyROS accepts models in which bounds are directly imposed on
Var
objects representing components of the variables \(z\)
and \(y\). These models are cast to
the form above
by reformulating the bounds as inequality constraints.
Note
A key requirement of PyROS is that each value of \(\left(x, z, q \right)\) maps to a unique value of \(y\), a property that is assumed to be properly enforced by the system of equality constraints \(\mathcal{J}\). If the mapping is not unique, then the selection of ‘state’ (i.e., not degree of freedom) variables \(y\) is incorrect, and one or more of the \(y\) variables should be appropriately redesignated to be part of either \(x\) or \(z\).
In order to cast the robust optimization counterpart of the deterministic model, we now assume that the uncertain parameters may attain any realization in a compact uncertainty set \(\mathcal{Q} \subseteq \mathbb{R}^{n_q}\) containing the nominal value \(q^{\text{nom}}\). The set \(\mathcal{Q}\) may be either continuous or discrete.
Based on the above notation, the form of the robust counterpart addressed by PyROS is
PyROS solves problems of this form using the Generalized Robust CuttingSet algorithm developed in [Isenberg_et_al].
When using PyROS, please consider citing the above paper.
PyROS Required Inputs
The required inputs to the PyROS solver are:
The deterministic optimization model
List of firststage (“design”) variables
List of secondstage (“control”) variables
List of parameters considered uncertain
The uncertainty set
Subordinate local and global nonlinear programming (NLP) solvers
These are more elaborately presented in the Solver Interface section.
Note
Any variables in the model not specified to be firststage or secondstage variables are automatically considered to be state variables.
PyROS Solver Interface
 class pyomo.contrib.pyros.PyROS[source]
PyROS (Pyomo Robust Optimization Solver) implementing a generalized robust cuttingset algorithm (GRCS) to solve twostage NLP optimization models under uncertainty.
 solve(model, first_stage_variables, second_stage_variables, uncertain_params, uncertainty_set, local_solver, global_solver, **kwds)[source]
Solve a model.
 Parameters:
model (ConcreteModel) – The deterministic model.
first_stage_variables (VarData, Var, or iterable of VarData/Var) – Firststage model variables (or design variables).
second_stage_variables (VarData, Var, or iterable of VarData/Var) – Secondstage model variables (or control variables).
uncertain_params (ParamData, Param, or iterable of ParamData/Param) – Uncertain model parameters. The mutable attribute for all uncertain parameter objects must be set to True.
uncertainty_set (UncertaintySet) – Uncertainty set against which the solution(s) returned will be confirmed to be robust.
local_solver (str or solver type) – Subordinate local NLP solver. If a str is passed, then the str is cast to
SolverFactory(local_solver)
.global_solver (str or solver type) – Subordinate global NLP solver. If a str is passed, then the str is cast to
SolverFactory(global_solver)
.
 Returns:
return_soln – Summary of PyROS termination outcome.
 Return type:
ROSolveResults
 Keyword Arguments:
time_limit (NonNegativeFloat, optional) – Wall time limit for the execution of the PyROS solver in seconds (including time spent by subsolvers). If None is provided, then no time limit is enforced.
keepfiles (bool, default=False) – Export subproblems with a nonacceptable termination status for debugging purposes. If True is provided, then the argument subproblem_file_directory must also be specified.
tee (bool, default=False) – Output subordinate solver logs for all subproblems.
load_solution (bool, default=True) – Load final solution(s) found by PyROS to the deterministic model provided.
symbolic_solver_labels (bool, default=False) – True to ensure the component names given to the subordinate solvers for every subproblem reflect the names of the corresponding Pyomo modeling components, False otherwise.
objective_focus (InEnum[ObjectiveType], default=<ObjectiveType.nominal: 2>) –
Objective focus for the master problems:
ObjectiveType.nominal: Optimize the objective function subject to the nominal uncertain parameter realization.
ObjectiveType.worst_case: Optimize the objective function subject to the worstcase uncertain parameter realization.
By default, ObjectiveType.nominal is chosen.
A worstcase objective focus is required for certification of robust optimality of the final solution(s) returned by PyROS. If a nominal objective focus is chosen, then only robust feasibility is guaranteed.
nominal_uncertain_param_vals (list, default=[]) – Nominal uncertain parameter realization. Entries should be provided in an order consistent with the entries of the argument uncertain_params. If an empty list is provided, then the values of the Param objects specified through uncertain_params are chosen.
decision_rule_order (In[0, 1, 2], default=0) –
Order (or degree) of the polynomial decision rule functions for approximating the adjustability of the second stage variables with respect to the uncertain parameters.
Choices are:
0: static recourse
1: affine recourse
2: quadratic recourse
solve_master_globally (bool, default=False) – True to solve all master problems with the subordinate global solver, False to solve all master problems with the subordinate local solver. Along with a worstcase objective focus (see argument objective_focus), solving the master problems to global optimality is required for certification of robust optimality of the final solution(s) returned by PyROS. Otherwise, only robust feasibility is guaranteed.
max_iter (positive int or 1, default=1) – Iteration limit. If 1 is provided, then no iteration limit is enforced.
robust_feasibility_tolerance (NonNegativeFloat, default=0.0001) – Relative tolerance for assessing maximal inequality constraint violations during the GRCS separation step.
separation_priority_order (dict, default={}) – Mapping from model inequality constraint names to positive integers specifying the priorities of their corresponding separation subproblems. A higher integer value indicates a higher priority. Constraints not referenced in the dict assume a priority of 0. Separation subproblems are solved in order of decreasing priority.
progress_logger (None, str or logging.Logger, default=<PreformattedLogger pyomo.contrib.pyros (INFO)>) – Logger (or name thereof) used for reporting PyROS solver progress. If None or a str is provided, then
progress_logger
is cast tologging.getLogger(progress_logger)
. In the default case, progress_logger is set to apyomo.contrib.pyros.util.PreformattedLogger
object of levellogging.INFO
.backup_local_solvers (str, solver type, or Iterable of str/solver type, default=[]) – Additional subordinate local NLP optimizers to invoke in the event the primary local NLP optimizer fails to solve a subproblem to an acceptable termination condition.
backup_global_solvers (str, solver type, or Iterable of str/solver type, default=[]) – Additional subordinate global NLP optimizers to invoke in the event the primary global NLP optimizer fails to solve a subproblem to an acceptable termination condition.
subproblem_file_directory (Path, optional) – Directory to which to export subproblems not successfully solved to an acceptable termination condition. In the event
keepfiles=True
is specified, a str or pathlike referring to an existing directory must be provided.bypass_local_separation (bool, default=False) – This is an advanced option. Solve all separation subproblems with the subordinate global solver(s) only. This option is useful for expediting PyROS in the event that the subordinate global optimizer(s) provided can quickly solve separation subproblems to global optimality.
bypass_global_separation (bool, default=False) – This is an advanced option. Solve all separation subproblems with the subordinate local solver(s) only. If True is chosen, then robustness of the final solution(s) returned by PyROS is not guaranteed, and a warning will be issued at termination. This option is useful for expediting PyROS in the event that the subordinate global optimizer provided cannot tractably solve separation subproblems to global optimality.
Note
Upon successful convergence of PyROS, the solution returned is certified to be robust optimal only if:
master problems are solved to global optimality (by specifying
solve_master_globally=True
)a worstcase objective focus is chosen (by specifying
objective_focus=ObjectiveType.worst_case
)
Otherwise, the solution returned is certified to only be robust feasible.
PyROS Uncertainty Sets
Uncertainty sets are represented by subclasses of
the UncertaintySet
abstract base class.
PyROS provides a suite of preimplemented subclasses representing
commonly used uncertainty sets.
Custom userdefined uncertainty set types may be implemented by
subclassing the
UncertaintySet
class.
The intersection of a sequence of concrete
UncertaintySet
instances can be easily constructed by instantiating the preimplemented
IntersectionSet
subclass.
The table that follows provides mathematical definitions of
the various abstract and preimplemented
UncertaintySet
subclasses.
Uncertainty Set Type 
Input Data 
Mathematical Definition 

\(\begin{array}{l} q ^{\text{L}} \in \mathbb{R}^{n}, \\ q^{\text{U}} \in \mathbb{R}^{n} \end{array}\) 
\(\{q \in \mathbb{R}^n \mid q^\mathrm{L} \leq q \leq q^\mathrm{U}\}\) 

\(\begin{array}{l} q^{0} \in \mathbb{R}^{n}, \\ \hat{q} \in \mathbb{R}_{+}^{n}, \\ \Gamma \in [0, n] \end{array}\) 
\(\left\{ q \in \mathbb{R}^{n} \middle \begin{array}{l} q = q^{0} + \hat{q} \circ \xi \\ \displaystyle \sum_{i=1}^{n} \xi_{i} \leq \Gamma \\ \xi \in [0, 1]^{n} \end{array} \right\}\) 

\(\begin{array}{l} q^{0} \in \mathbb{R}^{n}, \\ b \in \mathbb{R}_{+}^{L}, \\ B \in \{0, 1\}^{L \times n} \end{array}\) 
\(\left\{ q \in \mathbb{R}^{n} \middle \begin{array}{l} \begin{pmatrix} B \\ I \end{pmatrix} q \leq \begin{pmatrix} b + Bq^{0} \\ q^{0} \end{pmatrix} \end{array} \right\}\) 

\(\begin{array}{l} q^{0} \in \mathbb{R}^{n}, \\ \Psi \in \mathbb{R}^{n \times F}, \\ \beta \in [0, 1] \end{array}\) 
\(\left\{ q \in \mathbb{R}^{n} \middle \begin{array}{l} q = q^{0} + \Psi \xi \\ \displaystyle\bigg \sum_{j=1}^{F} \xi_{j} \bigg \leq \beta F \\ \xi \in [1, 1]^{F} \\ \end{array} \right\}\) 

\(\begin{array}{l} A \in \mathbb{R}^{m \times n}, \\ b \in \mathbb{R}^{m}\end{array}\) 
\(\{q \in \mathbb{R}^{n} \mid A q \leq b\}\) 

\(\begin{array}{l} q^0 \in \mathbb{R}^{n}, \\ \alpha \in \mathbb{R}_{+}^{n} \end{array}\) 
\(\left\{ q \in \mathbb{R}^{n} \middle \begin{array}{l} \displaystyle\sum_{\substack{i = 1: \\ \alpha_{i} > 0}}^{n} \left(\frac{q_{i}  q_{i}^{0}}{\alpha_{i}}\right)^2 \leq 1 \\ q_{i} = q_{i}^{0} \,\forall\,i : \alpha_{i} = 0 \end{array} \right\}\) 

\(\begin{array}{l} q^0 \in \mathbb{R}^n, \\ P \in \mathbb{S}_{++}^{n}, \\ s \in \mathbb{R}_{+} \end{array}\) 
\(\{q \in \mathbb{R}^{n} \mid (q  q^{0})^{\intercal} P^{1} (q  q^{0}) \leq s\}\) 

\(g: \mathbb{R}^{n} \to \mathbb{R}^{m}\) 
\(\{q \in \mathbb{R}^{n} \mid g(q) \leq 0\}\) 

\(q^{1}, q^{2},\dots , q^{S} \in \mathbb{R}^{n}\) 
\(\{q^{1}, q^{2}, \dots , q^{S}\}\) 

\(\mathcal{Q}_{1}, \mathcal{Q}_{2}, \dots , \mathcal{Q}_{m} \subset \mathbb{R}^{n}\) 
\(\displaystyle \bigcap_{i=1}^{m} \mathcal{Q}_{i}\) 
Note
Each of the PyROS uncertainty set classes inherits from the
UncertaintySet
abstract base class.
PyROS Uncertainty Set Classes
 class pyomo.contrib.pyros.uncertainty_sets.BoxSet(bounds)[source]
Bases:
UncertaintySet
A hyperrectangle (a.k.a. “box”).
 Parameters:
bounds ((N, 2) array_like) – Lower and upper bounds for each dimension of the set.
Examples
1D box set (interval):
>>> from pyomo.contrib.pyros import BoxSet >>> interval = BoxSet(bounds=[(1, 2)]) >>> interval.bounds array([[1, 2]])
2D box set:
>>> box_set = BoxSet(bounds=[[1, 2], [3, 4]]) >>> box_set.bounds array([[1, 2], [3, 4]])
5D hypercube with bounds 0 and 1 in each dimension:
>>> hypercube_5d = BoxSet(bounds=[[0, 1] for idx in range(5)]) >>> hypercube_5d.bounds array([[0, 1], [0, 1], [0, 1], [0, 1], [0, 1]])
 property bounds
Lower and upper bounds for each dimension of the set.
The bounds of a BoxSet instance can be changed, such that the dimension of the set remains unchanged.
 Type:
(N, 2) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the box set. This is numerically equivalent to the bounds attribute.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 class pyomo.contrib.pyros.uncertainty_sets.CardinalitySet(origin, positive_deviation, gamma)[source]
Bases:
UncertaintySet
A cardinalityconstrained (a.k.a. “gamma”) set.
 Parameters:
origin ((N,) array_like) – Origin of the set (e.g., nominal uncertain parameter values).
positive_deviation ((N,) array_like) – Maximal nonnegative coordinate deviation from the origin in each dimension.
gamma (numeric type) – Upper bound for the number of uncertain parameters which may realize their maximal deviations from the origin simultaneously.
Examples
A 3D cardinality set:
>>> from pyomo.contrib.pyros import CardinalitySet >>> gamma_set = CardinalitySet( ... origin=[0, 0, 0], ... positive_deviation=[1.0, 2.0, 1.5], ... gamma=1, ... ) >>> gamma_set.origin array([0, 0, 0]) >>> gamma_set.positive_deviation array([1. , 2. , 1.5]) >>> gamma_set.gamma 1
 property gamma
Upper bound for the number of uncertain parameters which may maximally deviate from their respective origin values simultaneously. Must be a numerical value ranging from 0 to the set dimension N.
Note that, mathematically, setting gamma to 0 reduces the set to a singleton containing the center, while setting gamma to the set dimension N makes the set mathematically equivalent to a BoxSet with bounds
numpy.array([origin, origin + positive_deviation]).T
. Type:
numeric type
 property origin
Origin of the cardinality set (e.g. nominal parameter values).
 Type:
(N,) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the cardinality set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)[source]
Determine whether a given point lies in the cardinality set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
True if the point lies in the set, False otherwise.
 Return type:
 property positive_deviation
Maximal coordinate deviations from the origin in each dimension. All entries are nonnegative.
 Type:
(N,) numpy.ndarray
 class pyomo.contrib.pyros.uncertainty_sets.BudgetSet(budget_membership_mat, rhs_vec, origin=None)[source]
Bases:
UncertaintySet
A budget set.
 Parameters:
budget_membership_mat ((L, N) array_like) – Incidence matrix of the budget constraints. Each row corresponds to a single budget constraint, and defines which uncertain parameters (which dimensions) participate in that row’s constraint.
rhs_vec ((L,) array_like) – Budget limits (upper bounds) with respect to the origin of the set.
origin ((N,) array_like or None, optional) – Origin of the budget set. If None is provided, then the origin is set to the zero vector.
Examples
3D budget set with one budget constraint and no origin chosen (hence origin defaults to 3D zero vector):
>>> from pyomo.contrib.pyros import BudgetSet >>> budget_set = BudgetSet( ... budget_membership_mat=[[1, 1, 1]], ... rhs_vec=[2], ... ) >>> budget_set.budget_membership_mat array([[1, 1, 1]]) >>> budget_set.budget_rhs_vec array([2]) >>> budget_set.origin array([0., 0., 0.])
3D budget set with two budget constraints and custom origin:
>>> budget_custom = BudgetSet( ... budget_membership_mat=[[1, 0, 1], [0, 1, 0]], ... rhs_vec=[1, 1], ... origin=[2, 2, 2], ... ) >>> budget_custom.budget_membership_mat array([[1, 0, 1], [0, 1, 0]]) >>> budget_custom.budget_rhs_vec array([1, 1]) >>> budget_custom.origin array([2, 2, 2])
 property budget_membership_mat
Incidence matrix of the budget constraints. Each row corresponds to a single budget constraint and defines which uncertain parameters participate in that row’s constraint.
 Type:
(L, N) numpy.ndarray
 property budget_rhs_vec
Budget limits (upper bounds) with respect to the origin.
 Type:
(L,) numpy.ndarray
 property coefficients_mat
Coefficient matrix of all polyhedral constraints defining the budget set. Composed from the incidence matrix used for defining the budget constraints and a coefficient matrix for individual uncertain parameter nonnegativity constraints.
This attribute cannot be set. The budget constraint incidence matrix may be altered through the budget_membership_mat attribute.
 Type:
(L + N, N) numpy.ndarray
 property origin
Origin of the budget set.
 Type:
(N,) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the budget set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 property rhs_vec
Righthand side vector for polyhedral constraints defining the budget set. This also includes entries for nonnegativity constraints on the uncertain parameters.
This attribute cannot be set, and is automatically determined given other attributes.
 Type:
(L + N,) numpy.ndarray
 class pyomo.contrib.pyros.uncertainty_sets.FactorModelSet(origin, number_of_factors, psi_mat, beta)[source]
Bases:
UncertaintySet
A factor model (a.k.a. “netalpha” model) set.
 Parameters:
origin ((N,) array_like) – Uncertain parameter values around which deviations are restrained.
number_of_factors (int) – Natural number representing the dimensionality of the space to which the set projects.
psi_mat ((N, F) array_like) – Matrix designating each uncertain parameter’s contribution to each factor. Each row is associated with a separate uncertain parameter. Each column is associated with a separate factor. Number of columns F of psi_mat should be equal to number_of_factors.
beta (numeric type) – Real value between 0 and 1 specifying the fraction of the independent factors that can simultaneously attain their extreme values.
Examples
A 4D factor model set with a 2D factor space:
>>> from pyomo.contrib.pyros import FactorModelSet >>> import numpy as np >>> fset = FactorModelSet( ... origin=np.zeros(4), ... number_of_factors=2, ... psi_mat=np.full(shape=(4, 2), fill_value=0.1), ... beta=0.5, ... ) >>> fset.origin array([0., 0., 0., 0.]) >>> fset.number_of_factors 2 >>> fset.psi_mat array([[0.1, 0.1], [0.1, 0.1], [0.1, 0.1], [0.1, 0.1]]) >>> fset.beta 0.5
 property beta
Real number ranging from 0 to 1 representing the fraction of the independent factors that can simultaneously attain their extreme values.
Note that, mathematically, setting
beta = 0
will enforce that as many factors will be above 0 as there will be below 0 (i.e., “zeronetalpha” model). Ifbeta = 1
, then the set is numerically equivalent to a BoxSet with bounds[origin  psi @ np.ones(F), origin + psi @ np.ones(F)].T
. Type:
numeric type
 property number_of_factors
Natural number representing the dimensionality F of the space to which the set projects.
This attribute is immutable, and may only be set at object construction. Typically, the number of factors is significantly less than the set dimension, but no restriction to that end is imposed here.
 Type:
 property origin
Uncertain parameter values around which deviations are restrained.
 Type:
(N,) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the factor model set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)[source]
Determine whether a given point lies in the factor model set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
True if the point lies in the set, False otherwise.
 Return type:
 property psi_mat
Matrix designating each uncertain parameter’s contribution to each factor. Each row is associated with a separate uncertain parameter. Each column with a separate factor.
 Type:
(N, F) numpy.ndarray
 class pyomo.contrib.pyros.uncertainty_sets.PolyhedralSet(lhs_coefficients_mat, rhs_vec)[source]
Bases:
UncertaintySet
A bounded convex polyhedron or polytope.
 Parameters:
lhs_coefficients_mat ((M, N) array_like) – Lefthand side coefficients for the linear inequality constraints defining the polyhedral set.
rhs_vec ((M,) array_like) – Righthand side values for the linear inequality constraints defining the polyhedral set. Each entry is an upper bound for the quantity
lhs_coefficients_mat @ x
, where x is an (N,) array representing any point in the polyhedral set.
Examples
2D polyhedral set with 4 defining inequalities:
>>> from pyomo.contrib.pyros import PolyhedralSet >>> pset = PolyhedralSet( ... lhs_coefficients_mat=[[1, 0], [0, 1], [1, 1], [1, 0]], ... rhs_vec=[0, 0, 0, 1], ... ) >>> pset.coefficients_mat array([[1, 0], [ 0, 1], [1, 1], [ 1, 0]]) >>> pset.rhs_vec array([0, 0, 0, 1])
 property coefficients_mat
Coefficient matrix for the (linear) inequality constraints defining the polyhedral set.
In tandem with the rhs_vec attribute, this matrix should be such that the polyhedral set is nonempty and bounded. Such a check is performed only at instance construction.
 Type:
(M, N) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the polyhedral set.
Currently, an empty list is returned, as the bounds cannot, in general, be computed without access to an optimization solver.
 point_in_set(point)
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 property rhs_vec
Righthand side values (upper bounds) for the (linear) inequality constraints defining the polyhedral set.
 Type:
(M,) numpy.ndarray
 class pyomo.contrib.pyros.uncertainty_sets.AxisAlignedEllipsoidalSet(center, half_lengths)[source]
Bases:
UncertaintySet
An axisaligned ellipsoid.
 Parameters:
center ((N,) array_like) – Center of the ellipsoid.
half_lengths ((N,) array_like) – Semiaxis lengths of the ellipsoid.
Examples
3D origincentered unit hypersphere:
>>> from pyomo.contrib.pyros import AxisAlignedEllipsoidalSet >>> sphere = AxisAlignedEllipsoidalSet( ... center=[0, 0, 0], ... half_lengths=[1, 1, 1] ... ) >>> sphere.center array([0, 0, 0]) >>> sphere.half_lengths array([1, 1, 1])
 property center
Center of the ellipsoid.
 Type:
(N,) numpy.ndarray
 property half_lengths
Semiaxis lengths.
 Type:
(N,) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the axisaligned ellipsoidal set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 class pyomo.contrib.pyros.uncertainty_sets.EllipsoidalSet(center, shape_matrix, scale=1)[source]
Bases:
UncertaintySet
A general ellipsoid.
 Parameters:
center ((N,) arraylike) – Center of the ellipsoid.
shape_matrix ((N, N) arraylike) – A positive definite matrix characterizing the shape and orientation of the ellipsoid.
scale (numeric type, optional) – Square of the factor by which to scale the semiaxes of the ellipsoid (i.e. the eigenvectors of the shape matrix). The default is 1.
Examples
3D origincentered unit hypersphere:
>>> from pyomo.contrib.pyros import EllipsoidalSet >>> import numpy as np >>> hypersphere = EllipsoidalSet( ... center=[0, 0, 0], ... shape_matrix=np.eye(3), ... scale=1, ... ) >>> hypersphere.center array([0, 0, 0]) >>> hypersphere.shape_matrix array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]) >>> hypersphere.scale 1
A 2D ellipsoid with custom rotation and scaling:
>>> rotated_ellipsoid = EllipsoidalSet( ... center=[1, 1], ... shape_matrix=[[4, 2], [2, 4]], ... scale=0.5, ... ) >>> rotated_ellipsoid.center array([1, 1]) >>> rotated_ellipsoid.shape_matrix array([[4, 2], [2, 4]]) >>> rotated_ellipsoid.scale 0.5
 property center
Center of the ellipsoid.
 Type:
(N,) numpy.ndarray
 property parameter_bounds
Bounds in each dimension of the ellipsoidal set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 property scale
Square of the factor by which to scale the semiaxes of the ellipsoid (i.e. the eigenvectors of the shape matrix).
 Type:
numeric type
 property shape_matrix
A positive definite matrix characterizing the shape and orientation of the ellipsoid.
 Type:
(N, N) numpy.ndarray
 class pyomo.contrib.pyros.uncertainty_sets.UncertaintySet[source]
Bases:
object
An object representing an uncertainty set to be passed to the PyROS solver.
An UncertaintySet object should be viewed as merely a container for data needed to parameterize the set it represents, such that the object’s attributes do not reference the components of a Pyomo modeling object.
 abstract property dim
Dimension of the uncertainty set (number of uncertain parameters in a corresponding optimization model of interest).
 abstract property parameter_bounds
Bounds for the value of each uncertain parameter constrained by the set (i.e. bounds for each set dimension).
 point_in_set(point)[source]
Determine whether a given point lies in the uncertainty set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
is_in_set – True if the point lies in the uncertainty set, False otherwise.
 Return type:
Notes
This method is invoked at the outset of a PyROS solver call to determine whether a userspecified nominal parameter realization lies in the uncertainty set.
 class pyomo.contrib.pyros.uncertainty_sets.DiscreteScenarioSet(scenarios)[source]
Bases:
UncertaintySet
A discrete set of finitely many uncertain parameter realizations (or scenarios).
 Parameters:
scenarios ((M, N) array_like) – A sequence of M distinct uncertain parameter realizations.
Examples
2D set with three scenarios:
>>> from pyomo.contrib.pyros import DiscreteScenarioSet >>> discrete_set = DiscreteScenarioSet( ... scenarios=[[1, 1], [2, 1], [1, 2]], ... ) >>> discrete_set.scenarios [(1, 1), (2, 1), (1, 2)]
 property parameter_bounds
Bounds in each dimension of the discrete scenario set.
 Returns:
List, length N, of 2tuples. Each tuple specifies the bounds in its corresponding dimension.
 Return type:
list of tuples
 point_in_set(point)[source]
Determine whether a given point lies in the discrete scenario set.
 Parameters:
point ((N,) arraylike) – Point (parameter value) of interest.
 Returns:
True if the point lies in the set, False otherwise.
 Return type:
 class pyomo.contrib.pyros.uncertainty_sets.IntersectionSet(**unc_sets)[source]
Bases:
UncertaintySet
An intersection of a sequence of uncertainty sets, each of which is represented by an UncertaintySet object.
 Parameters:
**unc_sets (dict) – PyROS UncertaintySet objects of which to construct an intersection. At least two uncertainty sets must be provided. All sets must be of the same dimension.
Examples
Intersection of origincentered 2D box (square) and 2D hypersphere (circle):
>>> from pyomo.contrib.pyros import ( ... BoxSet, AxisAlignedEllipsoidalSet, IntersectionSet, ... ) >>> square = BoxSet(bounds=[[1.5, 1.5], [1.5, 1.5]]) >>> circle = AxisAlignedEllipsoidalSet( ... center=[0, 0], ... half_lengths=[2, 2], ... ) >>> # to construct intersection, pass sets as keyword arguments >>> intersection = IntersectionSet(set1=square, set2=circle) >>> intersection.all_sets UncertaintySetList([...])
 property all_sets
List of the uncertainty sets of which to take the intersection. Must be of minimum length 2.
This attribute may be set through any iterable of UncertaintySet objects, and exhibits similar behavior to a list.
 Type:
UncertaintySetList
 property parameter_bounds
Uncertain parameter value bounds for the intersection set.
Currently, an empty list, as the bounds cannot, in general, be computed without access to an optimization solver.
PyROS Usage Example
In this section, we illustrate the usage of PyROS with a modeling example. The deterministic problem of interest is called hydro (available here), a QCQP taken from the GAMS Model Library. We have converted the model to Pyomo format using the GAMS Convert tool.
The hydro model features 31 variables,
of which 13 are degrees of freedom and 18 are state variables.
Moreover, there are
6 linear inequality constraints,
12 linear equality constraints,
6 nonlinear (quadratic) equality constraints,
and a quadratic objective.
We have extended this model by converting one objective coefficient,
two constraint coefficients, and one constraint righthand side
into Param
objects so that they can be considered uncertain later on.
Note
Per our analysis, the hydro problem satisfies the requirement that each value of \(\left(x, z, q \right)\) maps to a unique value of \(y\), which, in accordance with our earlier note, indicates a proper partitioning of the model variables into (firststage and secondstage) degrees of freedom and state variables.
Step 0: Import Pyomo and the PyROS Module
In anticipation of using the PyROS solver and building the deterministic Pyomo model:
>>> # === Required import ===
>>> import pyomo.environ as pyo
>>> import pyomo.contrib.pyros as pyros
>>> # === Instantiate the PyROS solver object ===
>>> pyros_solver = pyo.SolverFactory("pyros")
Step 1: Define the Deterministic Problem
The deterministic Pyomo model for hydro is shown below.
Note
Primitive data (Python literals) that have been hardcoded within a
deterministic model cannot be later considered uncertain,
unless they are first converted to Param
objects within
the ConcreteModel
object.
Furthermore, any Param
object that is to be later considered
uncertain must have the property mutable=True
.
Note
In case modifying the mutable
property inside the deterministic
model object itself is not straightforward in your context,
you may consider adding the following statement after
import pyomo.environ as pyo
but before defining the model
object: pyo.Param.DefaultMutable = True
.
For all Param
objects declared after this statement,
the attribute mutable
is set to True
by default.
Hence, nonmutable Param
objects are now declared by
explicitly passing the argument mutable=False
to the
Param
constructor.
>>> # === 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=8e5*m.x7**2 + m.x13 == 0)
>>> m.c15 = pyo.Constraint(expr=8e5*m.x8**2 + m.x14 == 0)
>>> m.c16 = pyo.Constraint(expr=8e5*m.x9**2 + m.x15 == 0)
>>> m.c17 = pyo.Constraint(expr=8e5*m.x10**2 + m.x16 == 0)
>>> m.c18 = pyo.Constraint(expr=8e5*m.x11**2 + m.x17 == 0)
>>> m.c19 = pyo.Constraint(expr=8e5*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)
Step 2: Define the Uncertainty
First, we need to collect into a list those Param
objects of our model
that represent potentially uncertain parameters.
For the purposes of our example, we shall assume uncertainty in the model
parameters [m.p[0], m.p[1], m.p[2], m.p[3]]
, for which we can
conveniently utilize the object m.p
(itself an indexed Param
object).
>>> # === Specify which parameters are uncertain ===
>>> # We can pass IndexedParams this way to PyROS,
>>> # or as an expanded list per index
>>> uncertain_parameters = [m.p]
Note
Any Param
object that is to be considered uncertain by PyROS
must have the property mutable=True
.
PyROS will seek to identify solutions that remain feasible for any
realization of these parameters included in an uncertainty set.
To that end, we need to construct an
UncertaintySet
object.
In our example, let us utilize the
BoxSet
constructor to specify
an uncertainty set of simple hyperrectangular geometry.
For this, we will assume each parameter value is uncertain within a
percentage of its nominal value. Constructing this specific
UncertaintySet
object can be done as follows:
>>> # === 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)
Step 3: Solve with PyROS
PyROS requires the user to supply one local and one global NLP solver to use for solving subproblems. For convenience, we shall have PyROS invoke BARON as both the local and the global NLP solver:
>>> # === Designate local and global NLP solvers ===
>>> local_solver = pyo.SolverFactory('baron')
>>> global_solver = pyo.SolverFactory('baron')
Note
Additional NLP optimizers can be automatically used in the event the primary
subordinate local or global optimizer passed
to the PyROS solve()
method
does not successfully solve a subproblem to an appropriate termination
condition. These alternative solvers are provided through the optional
keyword arguments backup_local_solvers
and backup_global_solvers
.
The final step in solving a model with PyROS is to construct the
remaining required inputs, namely
first_stage_variables
and second_stage_variables
.
Below, we present two separate cases.
PyROS Termination Conditions
PyROS will return one of six termination conditions upon completion.
These termination conditions are defined through the
pyrosTerminationCondition
enumeration
and tabulated below.
Termination Condition 
Description 


The final solution is robust optimal 

The final solution is robust feasible 

The posed problem is robust infeasible 

Maximum number of GRCS iteration reached 

Maximum number of time reached 

Unacceptable return status(es) from a usersupplied subsolver 
A SingleStage Problem
If we choose to designate all variables as either design or state variables,
without any control variables (i.e., all degrees of freedom are firststage),
we can use PyROS to solve the singlestage problem as shown below.
In particular, let us instruct PyROS that variables
m.x1
through m.x6
, m.x19
through m.x24
, and m.x31
correspond to firststage degrees of freedom.
>>> # === Designate which variables correspond to firststage
>>> # and secondstage 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,
... objective_focus=pyros.ObjectiveType.worst_case,
... solve_master_globally=True,
... load_solution=False,
... )
==============================================================================
PyROS: The Pyomo Robust Optimization Solver...
...

Robust optimal solution identified.

...

All done. Exiting PyROS.
==============================================================================
>>> # === 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(f"Final objective value: {single_stage_final_objective}")
Final objective value: 48367380.0
>>> print(f"PyROS termination condition: {termination_condition}")
PyROS termination condition: pyrosTerminationCondition.robust_optimal
PyROS Results Object
The results object returned by PyROS allows you to query the following information from the solve call:
iterations
: total iterations of the algorithmtime
: total wallclock time (or elapsed time) in secondspyros_termination_condition
: the GRCS algorithm termination conditionfinal_objective_value
: the final objective function value.
The preceding code snippet demonstrates how to retrieve this information.
If we pass load_solution=True
(the default setting)
to the solve()
method,
then the solution at which PyROS terminates will be loaded to
the variables of the original deterministic model.
Note that in the preceding code snippet,
we set load_solution=False
to ensure the next set of runs shown here can
utilize the initial point loaded to the original deterministic model,
as the initial point may affect the performance of subsolvers.
Note
The reported final_objective_value
and final model variable values
depend on the selection of the option objective_focus
.
The final_objective_value
is the sum of firststage
and secondstage objective functions.
If objective_focus = ObjectiveType.nominal
,
secondstage objective and variables are evaluated at
the nominal realization of the uncertain parameters, \(q^{\text{nom}}\).
If objective_focus = ObjectiveType.worst_case
, secondstage objective
and variables are evaluated at the worstcase realization
of the uncertain parameters, \(q^{k^\ast}\)
where \(k^\ast = \mathrm{argmax}_{k \in \mathcal{K}}~f_2(x,z^k,y^k,q^k)\).
A TwoStage Problem
For this next set of runs, we will
assume that some of the previously designated firststage degrees of
freedom are in fact secondstage degrees of freedom.
PyROS handles secondstage degrees of freedom via the use of polynomial
decision rules, of which the degree is controlled through the
optional keyword argument decision_rule_order
to the PyROS
solve()
method.
In this example, we select affine decision rules by setting
decision_rule_order=1
:
>>> # === Define the variable partitioning
>>> first_stage_variables =[m.x5, m.x6, m.x19, m.x22, m.x23, m.x24, m.x31]
>>> second_stage_variables = [m.x1, m.x2, m.x3, m.x4, m.x20, m.x21]
>>> # The remaining variables are implicitly designated to be state variables
>>> # === Call PyROS to solve the robust optimization problem ===
>>> results_2 = 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,
... objective_focus=pyros.ObjectiveType.worst_case,
... solve_master_globally=True,
... decision_rule_order=1,
... )
==============================================================================
PyROS: The Pyomo Robust Optimization Solver...
...

Robust optimal solution identified.

...

All done. Exiting PyROS.
==============================================================================
>>> # === Compare final objective to the singlestage solution
>>> two_stage_final_objective = round(
... pyo.value(results_2.final_objective_value),
... 1,
... )
>>> percent_difference = 100 * (
... two_stage_final_objective  single_stage_final_objective
... ) / (single_stage_final_objective)
>>> print("Percent objective change relative to constant decision rules "
... f"objective: {percent_difference:.2f}")
Percent objective change relative to constant decision rules objective: 24...
For this example, we notice a ~25% decrease in the final objective value when switching from a static decision rule (no secondstage recourse) to an affine decision rule.
Specifying Arguments Indirectly Through options
Like other Pyomo solver interface methods,
solve()
provides support for specifying options indirectly by passing
a keyword argument options
, whose value must be a dict
mapping names of arguments to solve()
to their desired values.
For example, the solve()
statement in the
twostage problem snippet
could have been equivalently written as:
>>> results_2 = 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,
... options={
... "objective_focus": pyros.ObjectiveType.worst_case,
... "solve_master_globally": True,
... "decision_rule_order": 1,
... },
... )
==============================================================================
PyROS: The Pyomo Robust Optimization Solver...
...

Robust optimal solution identified.

...

All done. Exiting PyROS.
==============================================================================
In the event an argument is passed directly
by position or keyword, and indirectly through options
,
an appropriate warning is issued,
and the value passed directly takes precedence over the value
passed through options
.
The Price of Robustness
In conjunction with standard Python control flow tools,
PyROS facilitates a “price of robustness” analysis for a model of interest
through the evaluation and comparison of the robust optimal
objective function value across any appropriately constructed hierarchy
of uncertainty sets.
In this example, we consider a sequence of
box uncertainty sets centered on the nominal uncertain
parameter realization, such that each box is parameterized
by a real value specifying a relative box size.
To this end, we construct an iterable called relative_deviation_list
whose entries are float
values representing the relative sizes.
We then loop through relative_deviation_list
so that for each relative
size, the corresponding robust optimal objective value
can be evaluated by creating an appropriate
BoxSet
instance and invoking the PyROS solver:
>>> # This takes a long time to run and therefore is not a doctest
>>> # === An array of maximum relative deviations from the nominal uncertain
>>> # parameter values to utilize in constructing box sets
>>> relative_deviation_list = [0.00, 0.10, 0.20, 0.30, 0.40]
>>> # === Final robust optimal objectives
>>> robust_optimal_objectives = []
>>> for relative_deviation in relative_deviation_list:
... bounds = [
... (nominal_values[i]  relative_deviation*nominal_values[i],
... nominal_values[i] + relative_deviation*nominal_values[i])
... for i in range(4)
... ]
... box_uncertainty_set = pyros.BoxSet(bounds = bounds)
... results = 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,
... objective_focus=pyros.ObjectiveType.worst_case,
... solve_master_globally=True,
... decision_rule_order=1,
... )
... is_robust_optimal = (
... results.pyros_termination_condition
... == pyros.pyrosTerminationCondition.robust_optimal
... )
... if not is_robust_optimal:
... print(f"Instance for relative deviation: {relative_deviation} "
... "not solved to robust optimality.")
... robust_optimal_objectives.append("")
... else:
... robust_optimal_objectives.append(str(results.final_objective_value))
For this example, we obtain the following price of robustness results:
Uncertainty Set Size (+/) ^{o} 
Robust Optimal Objective 
% Increase ^{x} 

0.00 
35,837,659.18 
0.00 % 
0.10 
36,135,182.66 
0.83 % 
0.20 
36,437,979.81 
1.68 % 
0.30 
43,478,190.91 
21.32 % 
0.40 

\(\text{}\) 
Notice that PyROS was successfully able to determine the robust infeasibility of the problem under the largest uncertainty set.
^{o} Relative Deviation from Nominal Realization
^{x} Relative to Deterministic Optimal Objective
This example clearly illustrates the potential impact of the uncertainty set size on the robust optimal objective function value and demonstrates the ease of implementing a price of robustness study for a given optimization problem under uncertainty.
PyROS Solver Log Output
The PyROS solver log output is controlled through the optional
progress_logger
argument, itself cast to
a standard Python logger (logging.Logger
) object
at the outset of a solve()
call.
The level of detail of the solver log output
can be adjusted by adjusting the level of the
logger object; see the following table.
Note that by default, progress_logger
is cast to a logger of level
logging.INFO
.
We refer the reader to the official Python logging library documentation for customization of Python logger objects; for a basic tutorial, see the logging HOWTO.
Logging Level 
Output Messages 








An example of an output log produced through the default PyROS progress logger is shown in the snippet that follows. Observe that the log contains the following information:
Introductory information (lines 1–18). Includes the version number, author information, (UTC) time at which the solver was invoked, and, if available, information on the local Git branch and commit hash.
Summary of solver options (lines 19–38).
Preprocessing information (lines 39–41). Wall time required for preprocessing the deterministic model and associated components, i.e. standardizing model components and adding the decision rule variables and equations.
Model component statistics (lines 42–58). Breakdown of model component statistics. Includes components added by PyROS, such as the decision rule variables and equations.
Iteration log table (lines 59–69). Summary information on the problem iterates and subproblem outcomes. The constituent columns are defined in detail in the table following the snippet.
Termination message (lines 70–71). Very brief summary of the termination outcome.
Timing statistics (lines 72–88). Tabulated breakdown of the solver timing statistics, based on a
pyomo.common.timing.HierarchicalTimer
printout. The identifiers are as follows:main
: Total time elapsed by the solver.main.dr_polishing
: Total time elapsed by the subordinate solvers on polishing of the decision rules.main.global_separation
: Total time elapsed by the subordinate solvers on global separation subproblems.main.local_separation
: Total time elapsed by the subordinate solvers on local separation subproblems.main.master
: Total time elapsed by the subordinate solvers on the master problems.main.master_feasibility
: Total time elapsed by the subordinate solvers on the master feasibility problems.main.preprocessing
: Total preprocessing time.main.other
: Total overhead time.
Termination statistics (lines 89–94). Summary of statistics related to the iterate at which PyROS terminates.
Exit message (lines 95–96).
1==============================================================================
2PyROS: The Pyomo Robust Optimization Solver, v1.2.11.
3 Pyomo version: 6.7.2
4 Commit hash: unknown
5 Invoked at UTC 20240328T00:00:00.000000
6
7Developed by: Natalie M. Isenberg (1), Jason A. F. Sherman (1),
8 John D. Siirola (2), Chrysanthos E. Gounaris (1)
9(1) Carnegie Mellon University, Department of Chemical Engineering
10(2) Sandia National Laboratories, Center for Computing Research
11
12The developers gratefully acknowledge support from the U.S. Department
13of Energy's Institute for the Design of Advanced Energy Systems (IDAES).
14==============================================================================
15================================= DISCLAIMER =================================
16PyROS is still under development.
17Please provide feedback and/or report any issues by creating a ticket at
18https://github.com/Pyomo/pyomo/issues/new/choose
19==============================================================================
20Solver options:
21 time_limit=None
22 keepfiles=False
23 tee=False
24 load_solution=True
25 symbolic_solver_labels=False
26 objective_focus=<ObjectiveType.worst_case: 1>
27 nominal_uncertain_param_vals=[0.13248000000000001, 4.97, 4.97, 1800]
28 decision_rule_order=1
29 solve_master_globally=True
30 max_iter=1
31 robust_feasibility_tolerance=0.0001
32 separation_priority_order={}
33 progress_logger=<PreformattedLogger pyomo.contrib.pyros (INFO)>
34 backup_local_solvers=[]
35 backup_global_solvers=[]
36 subproblem_file_directory=None
37 bypass_local_separation=False
38 bypass_global_separation=False
39 p_robustness={}
40
41Preprocessing...
42Done preprocessing; required wall time of 0.175s.
43
44Model statistics:
45 Number of variables : 62
46 Epigraph variable : 1
47 Firststage variables : 7
48 Secondstage variables : 6
49 State variables : 18
50 Decision rule variables : 30
51 Number of uncertain parameters : 4
52 Number of constraints : 81
53 Equality constraints : 24
54 Coefficient matching constraints : 0
55 Decision rule equations : 6
56 All other equality constraints : 18
57 Inequality constraints : 57
58 Firststage inequalities (incl. certain var bounds) : 10
59 Performance constraints (incl. var bounds) : 47
60
61Itn Objective 1Stg Shift 2Stg Shift #CViol Max Viol Wall Time (s)
62
630 3.5838e+07   5 1.8832e+04 1.741
641 3.5838e+07 3.5184e15 3.9404e15 10 4.2516e+06 3.766
652 3.5993e+07 1.8105e01 7.1406e01 13 5.2004e+06 6.288
663 3.6285e+07 5.1968e01 7.7753e01 4 1.7892e+04 8.247
674 3.6285e+07 9.1166e13 1.9702e15 0 7.1157e10g 11.456
68
69Robust optimal solution identified.
70
71Timing breakdown:
72
73Identifier ncalls cumtime percall %
74
75main 1 11.457 11.457 100.0
76 
77 dr_polishing 4 0.682 0.171 6.0
78 global_separation 47 1.109 0.024 9.7
79 local_separation 235 5.810 0.025 50.7
80 master 5 1.353 0.271 11.8
81 master_feasibility 4 0.247 0.062 2.2
82 preprocessing 1 0.429 0.429 3.7
83 other n/a 1.828 n/a 16.0
84 ======================================================
85===========================================================
86
87
88Termination stats:
89 Iterations : 5
90 Solve time (wall s) : 11.457
91 Final objective value : 3.6285e+07
92 Termination condition : pyrosTerminationCondition.robust_optimal
93
94All done. Exiting PyROS.
95==============================================================================
The iteration log table is designed to provide, in a concise manner, important information about the progress of the iterative algorithm for the problem of interest. The constituent columns are defined in the table that follows.
Column Name 
Definition 

Itn 
Iteration number. 
Objective 
Master solution objective function value.
If the objective of the deterministic model provided
has a maximization sense,
then the negative of the objective function value is displayed.
Expect this value to trend upward as the iteration number
increases.
If the master problems are solved globally
(by passing 
1Stg Shift 
Infinity norm of the relative difference between the firststage variable vectors of the master solutions of the current and previous iterations. Expect this value to trend downward as the iteration number increases. A dash (“”) is produced in lieu of a value if the current iteration number is 0, there are no firststage variables, or the master problem of the current iteration is not solved successfully. 
2Stg Shift 
Infinity norm of the relative difference between the secondstage variable vectors (evaluated subject to the nominal uncertain parameter realization) of the master solutions of the current and previous iterations. Expect this value to trend downward as the iteration number increases. A dash (“”) is produced in lieu of a value if the current iteration number is 0, there are no secondstage variables, or the master problem of the current iteration is not solved successfully. 
#CViol 
Number of performance constraints found to be violated during
the separation step of the current iteration.
Unless a custom prioritization of the model’s performance constraints
is specified (through the 
Max Viol 
Maximum scaled performance constraint violation. Expect this value to trend downward as the iteration number increases. A ‘g’ is appended to the value if the separation problems were solved globally during the current iteration. A dash (“”) is produced in lieu of a value if the separation routine is not invoked during the current iteration, or if there are no performance constraints. 
Wall time (s) 
Total time elapsed by the solver, in seconds, up to the end of the current iteration. 
Feedback and Reporting Issues
Please provide feedback and/or report any problems by opening an issue on the Pyomo GitHub page.