diff --git a/doc/OnlineDocs/explanation/solvers/pyros.rst b/doc/OnlineDocs/explanation/solvers/pyros.rst index bc5f59e9a1c..9311c8999fc 100644 --- a/doc/OnlineDocs/explanation/solvers/pyros.rst +++ b/doc/OnlineDocs/explanation/solvers/pyros.rst @@ -963,10 +963,10 @@ Observe that the log contains the following information: :linenos: ============================================================================== - PyROS: The Pyomo Robust Optimization Solver, v1.3.6. - Pyomo version: 6.9.2 - Commit hash: 41cd797e0 - Invoked at UTC 2025-03-13T16:20:31.105320+00:00 + PyROS: The Pyomo Robust Optimization Solver, v1.3.8. + Pyomo version: 6.9.3dev0 + Commit hash: unknown + Invoked at UTC 2025-05-05T00:00:00.000000+00:00 Developed by: Natalie M. Isenberg (1), Jason A. F. Sherman (1), John D. Siirola (2), Chrysanthos E. Gounaris (1) @@ -1025,11 +1025,10 @@ Observe that the log contains the following information: ------------------------------------------------------------------------------ Itn Objective 1-Stg Shift 2-Stg Shift #CViol Max Viol Wall Time (s) ------------------------------------------------------------------------------ - 0 3.5838e+07 - - 5 1.8832e+04 0.693 - 1 3.5838e+07 1.2289e-09 1.5876e-12 5 3.7762e+04 1.514 - 2 3.6129e+07 2.7244e-01 3.6878e-01 3 1.1093e+02 2.486 - 3 3.6269e+07 3.7352e-01 4.3227e-01 1 2.7711e+01 3.667 - 4 3.6285e+07 7.6526e-01 2.8426e-11 0 4.3364e-05g 6.291 + 0 3.5838e+07 - - 5 1.8832e+04 0.759 + 1 3.5838e+07 2.9329e-09 5.0030e-10 5 2.1295e+04 1.573 + 2 3.6285e+07 7.6526e-01 2.0398e-01 2 2.2457e+02 2.272 + 3 3.6285e+07 7.7212e-13 1.2525e-10 0 7.2940e-08g 5.280 ------------------------------------------------------------------------------ Robust optimal solution identified. ------------------------------------------------------------------------------ @@ -1037,22 +1036,22 @@ Observe that the log contains the following information: Identifier ncalls cumtime percall % ----------------------------------------------------------- - main 1 6.291 6.291 100.0 + main 1 5.281 5.281 100.0 ------------------------------------------------------ - dr_polishing 4 0.334 0.083 5.3 - global_separation 27 0.954 0.035 15.2 - local_separation 135 3.046 0.023 48.4 - master 5 1.027 0.205 16.3 - master_feasibility 4 0.133 0.033 2.1 - preprocessing 1 0.013 0.013 0.2 - other n/a 0.785 n/a 12.5 + dr_polishing 3 0.155 0.052 2.9 + global_separation 27 1.280 0.047 24.2 + local_separation 108 2.200 0.020 41.7 + master 4 0.727 0.182 13.8 + master_feasibility 3 0.103 0.034 1.9 + preprocessing 1 0.021 0.021 0.4 + other n/a 0.794 n/a 15.0 ====================================================== =========================================================== ------------------------------------------------------------------------------ Termination stats: - Iterations : 5 - Solve time (wall s) : 6.291 + Iterations : 4 + Solve time (wall s) : 5.281 Final objective value : 3.6285e+07 Termination condition : pyrosTerminationCondition.robust_optimal ------------------------------------------------------------------------------ @@ -1133,6 +1132,58 @@ The constituent columns are defined in the current iteration. +Separation Priority Ordering +---------------------------- +The PyROS solver supports custom prioritization of +the separation subproblems (and, thus, the constraints) +that are automatically derived from +a given model for robust optimization. +Users may specify separation priorities through: + +- (Recommended) :class:`~pyomo.core.base.suffix.Suffix` components + with local name ``pyros_separation_priority``, + declared on the model or any of its sub-blocks. + Each entry of every such + :class:`~pyomo.core.base.suffix.Suffix` + should map a + :class:`~pyomo.core.base.var.Var` + or :class:`~pyomo.core.base.constraint.Constraint` + component to a value that specifies the separation + priority of all constraints derived from that component +- The optional argument ``separation_priority_order`` + to the PyROS :py:meth:`~pyomo.contrib.pyros.pyros.PyROS.solve` + method. The argument should be castable to a :py:obj:`dict`, + of which each entry maps the full name of a + :class:`~pyomo.core.base.var.Var` + or :class:`~pyomo.core.base.constraint.Constraint` + component to a value that specifies the + separation priority of all constraints + derived from that component + +Specification via :class:`~pyomo.core.base.suffix.Suffix` components +takes precedence over specification via the solver argument +``separation_priority_order``. +Moreover, the precedence ordering among +:class:`~pyomo.core.base.suffix.Suffix` +components is handled by the Pyomo +:class:`~pyomo.core.base.suffix.SuffixFinder` utility. + +A separation priority can be either +a (real) number (i.e., of type :py:class:`int`, :py:class:`float`, etc.) +or :py:obj:`None`. +A higher number indicates a higher priority. +The default priority for all constraints is 0. +Therefore a constraint can be prioritized [or deprioritized] +over the default by mapping the constraint to a positive [or negative] number. +In practice, critical or dominant constraints are often +prioritized over algorithmic or implied constraints. + +Constraints that have been assigned a priority of :py:obj:`None` +are enforced subject to only the nominal uncertain parameter realization +provided by the user. Therefore, these constraints are not imposed robustly +and, in particular, are excluded from the separation problems. + + Feedback and Reporting Issues ------------------------------- Please provide feedback and/or report any problems by opening an issue on diff --git a/pyomo/contrib/pyros/CHANGELOG.txt b/pyomo/contrib/pyros/CHANGELOG.txt index e58b8648ef7..302d188c9e2 100644 --- a/pyomo/contrib/pyros/CHANGELOG.txt +++ b/pyomo/contrib/pyros/CHANGELOG.txt @@ -3,11 +3,24 @@ PyROS CHANGELOG =============== +------------------------------------------------------------------------------- +PyROS 1.3.8 28 Apr 2025 +------------------------------------------------------------------------------- +- Add Suffix-based interface for prioritizing separation problems +- Allow user to, through the separation priority ordering interface, + specify constraints that should be imposed subject to + the nominal realization only (and, therefore, not separated) +- Deprecate the optional argument `separation_priority_order`. + + ------------------------------------------------------------------------------- PyROS 1.3.7 06 Mar 2025 ------------------------------------------------------------------------------- - Modify reformulation of state-variable independent second-stage equality constraints for problems with discrete uncertainty sets +- Lift the iteration-wise DR efficiency for problems with DR variable-dependent + first-stage equality constraints + (such as those derived from coefficient matching) ------------------------------------------------------------------------------- diff --git a/pyomo/contrib/pyros/config.py b/pyomo/contrib/pyros/config.py index fb1e2001e8b..fde4f35375b 100644 --- a/pyomo/contrib/pyros/config.py +++ b/pyomo/contrib/pyros/config.py @@ -13,6 +13,7 @@ InEnum, Path, ) +from pyomo.common.deprecation import deprecation_warning from pyomo.common.errors import ApplicationError, PyomoException from pyomo.core.base import Var, VarData from pyomo.core.base.param import Param, ParamData @@ -58,6 +59,38 @@ def positive_int_or_minus_one(obj): positive_int_or_minus_one.domain_name = "positive int or -1" +def _deprecated_separation_priority_order(obj): + """ + Domain validator for argument `separation_priority_order`. + + As this argument has been deprecated, a deprecation warning + is issued through a WARNING-level logger message if + the argument is cast to a nonempty dict. + + Parameters + ---------- + obj : object + Argument value. + + Returns + ------- + separation_priority_order : dict + Argument value, cast to a dict. + """ + separation_priority_order = dict(obj) + if separation_priority_order: + deprecation_warning( + "The argument 'separation_priority_order' is deprecated. " + "Consider specifying separation priorities by declaring, on your " + "model, Suffix components with local name `pyros_separation_priority`.", + version="6.9.3dev0", + ) + return separation_priority_order + + +_deprecated_separation_priority_order.domain_name = dict.__name__ + + def uncertain_param_validator(uncertain_obj): """ Check that a component object modeling an @@ -720,17 +753,34 @@ def pyros_config(): "separation_priority_order", ConfigValue( default={}, - domain=dict, + domain=_deprecated_separation_priority_order, doc=( """ - 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. + (DEPRECATED) + A dict-like object, each entry of which + maps the full name of a model ``Var`` or ``Constraint`` + component to a value specifying the separation priority + for all constraints derived from the component. + A separation priority can be a numeric value or None. + A higher numeric value indicates a higher priority. + For all constraints, the default priority is 0. + (Inequality and equality) constraints with a + priority of None are excluded from + the separation problems and enforced subject to only + the nominal uncertain parameter realization in the master + problems. + Separation problems corresponding to inequality + constraints with numeric priorities are grouped by + priority. In every iteration, the groups are traversed + in descending order of priority, + until, within a group, constraint violations + are detected. + + *Deprecated since Pyomo 6.9.3dev0*: The argument + `separation_priority_order` is deprecated. + Specify separation priorities by declaring, on your + model, `Suffix` components with local name + 'pyros_separation_priority'. """ ), ), diff --git a/pyomo/contrib/pyros/pyros.py b/pyomo/contrib/pyros/pyros.py index ed4247ff89a..4da7b6672cd 100644 --- a/pyomo/contrib/pyros/pyros.py +++ b/pyomo/contrib/pyros/pyros.py @@ -33,7 +33,7 @@ ) -__version__ = "1.3.7" +__version__ = "1.3.8" default_pyros_solver_logger = setup_pyros_logger() diff --git a/pyomo/contrib/pyros/pyros_algorithm_methods.py b/pyomo/contrib/pyros/pyros_algorithm_methods.py index f944a310052..771411002e8 100644 --- a/pyomo/contrib/pyros/pyros_algorithm_methods.py +++ b/pyomo/contrib/pyros/pyros_algorithm_methods.py @@ -327,6 +327,7 @@ def ROSolver_iterative_solve(model_data): from_block=nominal_master_blk, clone_first_stage_components=False, ) + separation_data.points_added_to_master[(k + 1, 0)] = ( separation_results.violating_param_realization ) diff --git a/pyomo/contrib/pyros/separation_problem_methods.py b/pyomo/contrib/pyros/separation_problem_methods.py index 5d1953a5a51..be928c83d18 100644 --- a/pyomo/contrib/pyros/separation_problem_methods.py +++ b/pyomo/contrib/pyros/separation_problem_methods.py @@ -387,22 +387,44 @@ def group_ss_ineq_constraints_by_priority(separation_data): Keys are sorted in descending order (i.e. highest priority first). """ + separation_data.config.progress_logger.debug( + "Grouping second-stage inequality constraints by separation priority..." + ) + ss_ineq_cons = separation_data.separation_model.second_stage.inequality_cons separation_priority_groups = dict() for name, ss_ineq_con in ss_ineq_cons.items(): - # by default, priority set to 0 priority = separation_data.separation_priority_order[name] cons_with_same_priority = separation_priority_groups.setdefault(priority, []) cons_with_same_priority.append(ss_ineq_con) # sort separation priority groups - return { + numeric_priority_grp_items = [ + (priority, cons) for priority, cons in separation_priority_groups.items() + ] + sorted_priority_groups = { priority: ss_ineq_cons - for priority, ss_ineq_cons in sorted( - separation_priority_groups.items(), reverse=True - ) + for priority, ss_ineq_cons in sorted(numeric_priority_grp_items, reverse=True) } + num_priority_groups = len(sorted_priority_groups) + separation_data.config.progress_logger.debug( + f"Found {num_priority_groups} separation " + f"priority group{'s' if num_priority_groups != 1 else ''}." + ) + separation_data.config.progress_logger.debug( + "Separation priority grouping statistics:" + ) + separation_data.config.progress_logger.debug( + f" {'Priority':20s}{'# Ineq Cons':15s}" + ) + for priority, cons in sorted_priority_groups.items(): + separation_data.config.progress_logger.debug( + f" {priority:<20d}{len(cons):<15d}" + ) + + return sorted_priority_groups + def get_worst_discrete_separation_solution( ss_ineq_con, config, ss_ineq_cons_to_evaluate, discrete_solve_results @@ -567,7 +589,7 @@ def perform_separation_loop(separation_data, master_data, solve_globally): master_data=master_data, ss_ineq_cons=all_ss_ineq_constraints, ) - sorted_priority_groups = group_ss_ineq_constraints_by_priority(separation_data) + sorted_priority_groups = separation_data.separation_priority_groups uncertainty_set_is_discrete = ( config.uncertainty_set.geometry == Geometry.DISCRETE_SCENARIOS ) @@ -629,11 +651,12 @@ def perform_separation_loop(separation_data, master_data, solve_globally): all_solve_call_results = ComponentMap() priority_groups_enum = enumerate(sorted_priority_groups.items()) + solve_adverb = "Globally" if solve_globally else "Locally" for group_idx, (priority, ss_ineq_constraints) in priority_groups_enum: priority_group_solve_call_results = ComponentMap() + for idx, ss_ineq_con in enumerate(ss_ineq_constraints): # log progress of separation loop - solve_adverb = "Globally" if solve_globally else "Locally" config.progress_logger.debug( f"{solve_adverb} separating second-stage inequality constraint " f"{get_con_name_repr(separation_data.separation_model, ss_ineq_con)} " @@ -1294,6 +1317,8 @@ def __init__(self, model_data): else: self.idxs_of_master_scenarios = None + self.separation_priority_groups = group_ss_ineq_constraints_by_priority(self) + def solve_separation(self, master_data): """ Solve the separation problem. diff --git a/pyomo/contrib/pyros/tests/test_grcs.py b/pyomo/contrib/pyros/tests/test_grcs.py index e4ce53df14c..74cc7ad3039 100644 --- a/pyomo/contrib/pyros/tests/test_grcs.py +++ b/pyomo/contrib/pyros/tests/test_grcs.py @@ -33,7 +33,7 @@ ) from pyomo.common.errors import ApplicationError, InfeasibleConstraintException from pyomo.core.expr import replace_expressions -from pyomo.environ import maximize as pyo_max, units as u +from pyomo.environ import assert_optimal_termination, maximize as pyo_max, units as u from pyomo.opt import ( SolverResults, SolverStatus, @@ -51,6 +51,7 @@ Objective, Param, SolverFactory, + Suffix, Var, exp, log, @@ -2050,6 +2051,209 @@ def test_two_stage_set_nonstatic_dr_robust_opt(self, use_discrete_set, dr_order) self.assertAlmostEqual(m.z.value, 2, places=4) +@unittest.skipUnless(ipopt_available, "IPOPT not available.") +class TestPyROSSeparationPriorityOrder(unittest.TestCase): + """ + Test PyROS solver behavior with respect to specification + of separation priorities. + """ + + def test_priority_nominal_only_eq(self): + m = ConcreteModel() + m.q = Param(initialize=0, mutable=True) + m.x = Var(bounds=[-2, 2]) + m.z = Var(bounds=(None, m.q)) + m.eq_con = Constraint(expr=m.z == m.q**2) + m.obj = Objective(expr=m.x + m.z, sense=minimize) + m.pyros_separation_priority = Suffix() + # enforce equality only nominally, or else model would be + # robust infeasible with [0, 1] interval uncertainty set + # due to coefficient matching of the equality + m.pyros_separation_priority[m.eq_con] = None + pyros_solver = SolverFactory("pyros") + ipopt = SolverFactory("ipopt") + res = pyros_solver.solve( + model=m, + first_stage_variables=[m.x], + second_stage_variables=[m.z], + uncertain_params=[m.q], + uncertainty_set=BoxSet([[0, 1]]), + local_solver=ipopt, + global_solver=ipopt, + objective_focus="worst_case", + bypass_global_separation=True, + solve_master_globally=True, + decision_rule_order=0, + ) + self.assertEqual( + res.pyros_termination_condition, pyrosTerminationCondition.robust_optimal + ) + self.assertEqual(m.x.value, -2) + self.assertEqual(m.z.value, 0) + self.assertAlmostEqual(res.final_objective_value, -2, places=4) + # z is essentially fixed due to the equality, + # and x not involved in any constraints, so: + self.assertEqual(res.iterations, 1) + + def test_priority_nominal_only_var_bounds(self): + m = ConcreteModel() + m.q = Param(initialize=0, mutable=True) + m.x = Var(bounds=[-2, 2]) + m.y = Var(bounds=(m.q, None)) + m.eq_con = Constraint(expr=m.y == m.q**2) + m.obj = Objective(expr=m.x + m.y, sense=minimize) + m.pyros_separation_priority = Suffix() + # enforce bounds only nominally, or else model is robust + # infeasible with [0, 1] interval uncertainty set + m.pyros_separation_priority[m.y] = None + pyros_solver = SolverFactory("pyros") + ipopt = SolverFactory("ipopt") + res = pyros_solver.solve( + model=m, + first_stage_variables=[m.x], + second_stage_variables=[], + uncertain_params=[m.q], + uncertainty_set=BoxSet([[0, 1]]), + local_solver=ipopt, + global_solver=ipopt, + objective_focus="worst_case", + bypass_global_separation=True, + solve_master_globally=True, + ) + self.assertEqual( + res.pyros_termination_condition, pyrosTerminationCondition.robust_optimal + ) + self.assertEqual(m.x.value, -2) + self.assertEqual(m.y.value, 1) + # epigraph constraint is separated (due to worst-case focus), + # need only one more iteration to achieve robustness + self.assertEqual(res.iterations, 2) + + def test_priority_nominal_only_ineq(self): + m = ConcreteModel() + m.q = Param(initialize=0, mutable=True) + m.x = Var(bounds=[-2, 2]) + m.y = Var() + m.con = Constraint(expr=m.y >= m.q) + m.eq_con = Constraint(expr=m.y == m.q**2) + m.obj = Objective(expr=m.x + m.y, sense=minimize) + m.pyros_separation_priority = Suffix() + # enforce inequality only nominally, or else model is robust + # infeasible with [0, 1] interval uncertainty set + m.pyros_separation_priority[m.con] = None + pyros_solver = SolverFactory("pyros") + ipopt = SolverFactory("ipopt") + res = pyros_solver.solve( + model=m, + first_stage_variables=[m.x], + second_stage_variables=[], + uncertain_params=[m.q], + uncertainty_set=BoxSet([[0, 1]]), + local_solver=ipopt, + global_solver=ipopt, + objective_focus="worst_case", + bypass_global_separation=True, + solve_master_globally=True, + ) + self.assertEqual( + res.pyros_termination_condition, pyrosTerminationCondition.robust_optimal + ) + self.assertEqual(m.x.value, -2) + self.assertEqual(m.y.value, 1) + + def test_priority_skip_all_separation(self): + m = build_leyffer_two_cons() + m_det = m.clone() + m.pyros_separation_priority = Suffix() + m.pyros_separation_priority[None] = None + interval = BoxSet(bounds=[(0.25, 2)]) + pyros_solver = SolverFactory("pyros") + local_subsolver = SolverFactory('ipopt') + global_subsolver = SolverFactory("ipopt") + + res = pyros_solver.solve( + model=m, + first_stage_variables=[m.x1], + second_stage_variables=[m.x2], + uncertain_params=[m.u], + uncertainty_set=interval, + local_solver=local_subsolver, + global_solver=global_subsolver, + objective_focus="worst_case", + bypass_global_separation=True, + # note: this gets overridden by the priority suffix, + # and is therefore ignored + separation_priority_order={"con1": 2}, + decision_rule_order=1, + ) + + self.assertEqual( + res.pyros_termination_condition, + pyrosTerminationCondition.robust_feasible, + msg="Returned termination condition is not return robust_optimal.", + ) + self.assertEqual(res.iterations, 1) + assert_optimal_termination(local_subsolver.solve(m_det)) + # when all separation problems bypassed, PyROS reduces to a + # solving the deterministic model + self.assertAlmostEqual(m.x1.value, m_det.x1.value, places=4) + self.assertAlmostEqual(m.x2.value, m_det.x2.value, places=4) + self.assertAlmostEqual(m.x3.value, m_det.x3.value, places=4) + self.assertAlmostEqual(value(m.obj), value(m_det.obj), places=4) + self.assertAlmostEqual(res.final_objective_value, value(m_det.obj), places=4) + + def test_priority_order_invariant(self): + m = build_leyffer_two_cons() + m2 = m.clone() + interval = BoxSet(bounds=[(0.25, 2)]) + pyros_solver = SolverFactory("pyros") + local_subsolver = SolverFactory('ipopt') + global_subsolver = SolverFactory("ipopt") + res1 = pyros_solver.solve( + model=m, + first_stage_variables=[m.x1], + second_stage_variables=[m.x2], + uncertain_params=[m.u], + uncertainty_set=interval, + local_solver=local_subsolver, + global_solver=global_subsolver, + objective_focus="worst_case", + bypass_global_separation=True, + separation_priority_order={"con1": 2}, + ) + self.assertEqual( + res1.pyros_termination_condition, + pyrosTerminationCondition.robust_feasible, + msg="Returned termination condition is not return robust_optimal.", + ) + + m2.pyros_separation_priority = Suffix() + m2.pyros_separation_priority[m2.con1] = 2 + res2 = pyros_solver.solve( + model=m2, + first_stage_variables=[m2.x1], + second_stage_variables=[m2.x2], + uncertain_params=[m2.u], + uncertainty_set=interval, + local_solver=local_subsolver, + global_solver=global_subsolver, + objective_focus="worst_case", + bypass_global_separation=True, + ) + self.assertEqual( + res2.pyros_termination_condition, + pyrosTerminationCondition.robust_feasible, + msg="Returned termination condition is not return robust_optimal.", + ) + + # confirm results are identical + self.assertEqual(res2.iterations, res1.iterations) + self.assertEqual(res2.final_objective_value, res1.final_objective_value) + self.assertEqual(m.x1.value, m2.x1.value) + self.assertEqual(m.x2.value, m2.x2.value) + self.assertEqual(m.x3.value, m2.x3.value) + + @unittest.skipUnless(baron_available, "BARON not available") class TestReformulateSecondStageEqualitiesDiscrete(unittest.TestCase): """ diff --git a/pyomo/contrib/pyros/tests/test_preprocessor.py b/pyomo/contrib/pyros/tests/test_preprocessor.py index 8997e2a2707..ceccb14000d 100644 --- a/pyomo/contrib/pyros/tests/test_preprocessor.py +++ b/pyomo/contrib/pyros/tests/test_preprocessor.py @@ -34,6 +34,7 @@ RangeSet, maximize, Block, + Suffix, ) from pyomo.core.base.set_types import NonNegativeReals, NonPositiveReals, Reals from pyomo.core.expr import ( @@ -67,6 +68,7 @@ VariablePartitioning, preprocess_model_data, log_model_statistics, + DEFAULT_SEPARATION_PRIORITY, ) parameterized, param_available = attempt_import('parameterized') @@ -108,8 +110,11 @@ def build_simple_test_model_data(self): m.c4 = Constraint(expr=m.x2**2 + m.y[1] + m.y[2] + m.y[3] + m.y[4] == 0) m.c5 = Constraint(expr=m.x2 + 2 * m.y[2] + m.y[3] + 2 * m.y[4] == 0) - model_data = Bunch() - model_data.config = Bunch() + model_data = ModelData( + original_model=m, + config=Bunch(separation_priority_order=dict()), + timing=None, + ) model_data.working_model = ConcreteModel() model_data.working_model.user_model = mdl = m.clone() model_data.working_model.uncertain_params = [mdl.q, mdl.q2] @@ -334,9 +339,7 @@ def build_test_model_data(self): """ Build model data object for the preprocessor. """ - model_data = Bunch() - model_data.config = Bunch() - model_data.original_model = m = ConcreteModel() + m = ConcreteModel() # PARAMS: one uncertain, one certain m.p = Param(initialize=2, mutable=True) @@ -415,6 +418,8 @@ def build_test_model_data(self): state_variables=[m.y1, m.y2], ) + model_data = ModelData(original_model=m, config=Bunch(), timing=None) + return model_data, user_var_partitioning def test_setup_working_model(self): @@ -640,12 +645,7 @@ def build_simple_test_model_data(self): Build simple model data object for turning bounds to constraints. """ - model_data = Bunch() - model_data.config = Bunch() - - model_data.working_model = ConcreteModel() - model_data.working_model.user_model = m = ConcreteModel() - + m = ConcreteModel() m.q1 = Param(initialize=1, mutable=True) m.q2 = Param(initialize=1, mutable=True) m.p1 = Param(initialize=5, mutable=True) @@ -661,7 +661,16 @@ def build_simple_test_model_data(self): m.z8 = Var(domain=RangeSet(0, 5, 0), bounds=[m.q1, m.q2]) m.z9 = Var(domain=RangeSet(0, 5, 0), bounds=[m.q1, m.p1]) m.z10 = Var(domain=RangeSet(0, 5, 0), bounds=[m.q1, m.p2]) + m.z11 = Var(bounds=[m.q1, m.q1]) + model_data = ModelData( + original_model=None, + config=Bunch(separation_priority_order=dict()), + timing=None, + ) + + model_data.working_model = ConcreteModel() + model_data.working_model.user_model = m model_data.working_model.uncertain_params = [m.q1, m.q2, m.p1] model_data.working_model.effective_uncertain_params = [m.q1, m.q2] @@ -689,11 +698,32 @@ def test_turn_nonadjustable_bounds_to_constraints(self): # mock effective partitioning for testing ep = model_data.working_model.effective_var_partitioning = Bunch() - ep.first_stage_variables = [m.z1, m.z2, m.z3, m.z4, m.z5, m.z6, m.z7, m.z8] + ep.first_stage_variables = [ + m.z1, + m.z2, + m.z3, + m.z4, + m.z5, + m.z6, + m.z7, + m.z8, + m.z11, + ] ep.second_stage_variables = [m.z9] ep.state_variables = [m.z10] effective_first_stage_var_set = ComponentSet(ep.first_stage_variables) + # also want to test resolution of separation priorities + model_data.config.separation_priority_order["z3"] = 10 + model_data.config.separation_priority_order["z8"] = 9 + model_data.config.separation_priority_order["z11"] = None + m.pyros_separation_priority = Suffix() + m.pyros_separation_priority[m.z4] = 1 + m.pyros_separation_priority[m.z6] = 2 + # note: this suffix entry, rather than the + # config specification, should determine the priority + m.pyros_separation_priority[m.z8] = 4 + original_var_domains_and_bounds = ComponentMap( (var, (var.domain, get_var_bound_pairs(var)[1])) for var in model_data.working_model.user_model.component_data_objects(Var) @@ -710,6 +740,7 @@ def test_turn_nonadjustable_bounds_to_constraints(self): (m.z6, ((None, None), ["eq"])), (m.z7, ((None, None), ["eq"])), (m.z8, ((None, None), ["lower", "upper"])), + (m.z11, ((m.q1, m.q1), [])), ) ) @@ -821,15 +852,30 @@ def test_turn_nonadjustable_bounds_to_constraints(self): ) # check separation priorities - for con_name in second_stage.inequality_cons: - self.assertEqual( - model_data.separation_priority_order[con_name], - 0, - msg=( - f"Separation priority for entry {con_name!r} of second-stage " - "inequalities not as expected." - ), - ) + self.assertEqual( + len(model_data.separation_priority_order), + (len(second_stage.inequality_cons) + len(second_stage.equality_cons)), + ) + self.assertEqual( + model_data.separation_priority_order["var_z4_uncertain_lower_bound_con"], 1 + ) + self.assertEqual( + model_data.separation_priority_order["var_z5_uncertain_upper_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z6_uncertain_eq_bound_con"], 2 + ) + self.assertEqual( + model_data.separation_priority_order["var_z7_uncertain_eq_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_uncertain_lower_bound_con"], 4 + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_uncertain_upper_bound_con"], 4 + ) def test_turn_adjustable_bounds_to_constraints(self): """ @@ -845,13 +891,19 @@ def test_turn_adjustable_bounds_to_constraints(self): """ model_data = self.build_simple_test_model_data() + model_data.working_model.first_stage = Block() + model_data.working_model.first_stage.inequality_cons = Constraint(Any) + model_data.working_model.first_stage.dr_independent_equality_cons = Constraint( + Any + ) + m = model_data.working_model.user_model # simple mock partitioning for the test ep = model_data.working_model.effective_var_partitioning = Bunch() ep.first_stage_variables = [m.z9, m.z10] ep.second_stage_variables = [m.z1, m.z2, m.z3, m.z4, m.z5, m.z6] - ep.state_variables = [m.z7, m.z8] + ep.state_variables = [m.z7, m.z8, m.z11] effective_first_stage_var_set = ComponentSet(ep.first_stage_variables) original_var_domains_and_bounds = ComponentMap( @@ -859,6 +911,16 @@ def test_turn_adjustable_bounds_to_constraints(self): for var in model_data.working_model.user_model.component_data_objects(Var) ) + model_data.config.separation_priority_order["z3"] = 10 + model_data.config.separation_priority_order["z8"] = 9 + model_data.config.separation_priority_order["z11"] = None + m.pyros_separation_priority = Suffix() + m.pyros_separation_priority[m.z4] = 1 + m.pyros_separation_priority[m.z6] = 2 + # note: this suffix entry, rather than the + # config specification, should determine the priority + m.pyros_separation_priority[m.z8] = 4 + turn_adjustable_var_bounds_to_constraints(model_data) for var, (orig_domain, orig_bounds) in original_var_domains_and_bounds.items(): @@ -926,8 +988,11 @@ def test_turn_adjustable_bounds_to_constraints(self): ), ) + first_stage = model_data.working_model.first_stage second_stage = model_data.working_model.second_stage + self.assertEqual(len(first_stage.dr_independent_equality_cons), 1) + self.assertEqual(len(first_stage.inequality_cons), 0) self.assertEqual(len(second_stage.inequality_cons), 10) self.assertEqual(len(second_stage.equality_cons), 5) @@ -1007,17 +1072,66 @@ def test_turn_adjustable_bounds_to_constraints(self): second_stage.inequality_cons["var_z8_uncertain_upper_bound_con"].expr, m.z8 <= m.q2, ) + fs_dr_indep_cons = first_stage.dr_independent_equality_cons + assertExpressionsEqual( + self, fs_dr_indep_cons["var_z11_uncertain_eq_bound_con"].expr, m.z11 == m.q1 + ) # check separation priorities - for con_name in second_stage.inequality_cons: - self.assertEqual( - model_data.separation_priority_order[con_name], - 0, - msg=( - f"Separation priority for entry {con_name!r} of second-stage " - "inequalities not as expected." - ), - ) + self.assertEqual( + len(model_data.separation_priority_order), + (len(second_stage.inequality_cons) + len(second_stage.equality_cons)), + ) + self.assertEqual( + model_data.separation_priority_order["var_z2_certain_eq_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z3_certain_lower_bound_con"], 10 + ) + self.assertEqual( + model_data.separation_priority_order["var_z3_certain_upper_bound_con"], 10 + ) + self.assertEqual( + model_data.separation_priority_order["var_z4_certain_eq_bound_con"], 1 + ) + self.assertEqual( + model_data.separation_priority_order["var_z4_uncertain_lower_bound_con"], 1 + ) + self.assertEqual( + model_data.separation_priority_order["var_z5_certain_eq_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z5_uncertain_upper_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z6_certain_lower_bound_con"], 2 + ) + self.assertEqual( + model_data.separation_priority_order["var_z6_uncertain_eq_bound_con"], 2 + ) + self.assertEqual( + model_data.separation_priority_order["var_z7_certain_upper_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z7_uncertain_eq_bound_con"], + DEFAULT_SEPARATION_PRIORITY, + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_certain_lower_bound_con"], 4 + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_certain_upper_bound_con"], 4 + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_uncertain_lower_bound_con"], 4 + ) + self.assertEqual( + model_data.separation_priority_order["var_z8_uncertain_upper_bound_con"], 4 + ) class TestStandardizeInequalityConstraints(unittest.TestCase): @@ -1030,8 +1144,7 @@ def build_simple_test_model_data(self): Build model data object for testing constraint standardization routines. """ - model_data = Bunch() - model_data.config = Bunch() + model_data = ModelData(original_model=None, timing=None, config=Bunch()) model_data.working_model = ConcreteModel() model_data.working_model.user_model = m = Block() @@ -1058,6 +1171,7 @@ def build_simple_test_model_data(self): m.c10 = Constraint(expr=m.y1 <= m.q**2) m.c11 = Constraint(expr=m.z2 <= m.q) m.c12 = Constraint(expr=(m.q**2, m.x1, sin(m.p) * m.q_cert)) + m.c13 = Constraint(expr=m.x1 <= m.q) m.c11.deactivate() @@ -1082,6 +1196,7 @@ def build_simple_test_model_data(self): m.c9, m.c10, m.c12, + m.c13, ] ep = model_data.working_model.effective_var_partitioning = Bunch() @@ -1102,13 +1217,19 @@ def test_standardize_inequality_constraints(self): m = working_model.user_model model_data.config.separation_priority_order = dict(c3=1, c5=2) + m.pyros_separation_priority = Suffix() + m.pyros_separation_priority[m.c5] = 10 + m.pyros_separation_priority[m.c12] = 5 + m.pyros_separation_priority[m.c13] = None standardize_inequality_constraints(model_data) fs_ineq_cons = working_model.first_stage.inequality_cons ss_ineq_cons = working_model.second_stage.inequality_cons + sep_priority_dict = model_data.separation_priority_order - self.assertEqual(len(fs_ineq_cons), 4) + self.assertEqual(len(fs_ineq_cons), 5) self.assertEqual(len(ss_ineq_cons), 13) + self.assertEqual(len(sep_priority_dict), 13) self.assertFalse(m.c1.active) new_c1_con = fs_ineq_cons["ineq_con_c1"] @@ -1128,7 +1249,7 @@ def test_standardize_inequality_constraints(self): new_c3_con = ss_ineq_cons["ineq_con_c3_lower_bound_con"] self.assertTrue(new_c3_con.active) assertExpressionsEqual(self, new_c3_con.expr, -m.x1 <= -m.q) - self.assertEqual(model_data.separation_priority_order[new_c3_con.index()], 1) + self.assertEqual(sep_priority_dict[new_c3_con.index()], 1) # m.x1 - 2 * m.q <= 0; # single second-stage inequality. modify in place @@ -1158,12 +1279,8 @@ def test_standardize_inequality_constraints(self): self.assertTrue(new_c5_lower_bound_con.active) assertExpressionsEqual(self, new_c5_lower_bound_con.expr, -m.x2 <= -m.q) assertExpressionsEqual(self, new_c5_upper_bound_con.expr, m.x2 <= 2 * m.q) - self.assertEqual( - model_data.separation_priority_order[new_c5_lower_bound_con.index()], 2 - ) - self.assertEqual( - model_data.separation_priority_order[new_c5_upper_bound_con.index()], 2 - ) + self.assertEqual(sep_priority_dict[new_c5_lower_bound_con.index()], 10) + self.assertEqual(sep_priority_dict[new_c5_upper_bound_con.index()], 10) # single second-stage inequality self.assertFalse(m.c6.active) @@ -1219,13 +1336,23 @@ def test_standardize_inequality_constraints(self): assertExpressionsEqual( self, new_c12_upper_bound_con.expr, m.x1 <= sin(m.p) * m.q_cert ) + self.assertEqual(sep_priority_dict[new_c12_lower_bound_con.index()], 5) + + self.assertFalse(m.c13.active) + new_c13_upper_bound_con = fs_ineq_cons["ineq_con_c13"] + self.assertTrue(new_c13_upper_bound_con.active) + assertExpressionsEqual(self, new_c13_upper_bound_con.expr, m.x1 <= m.q) # check separation priorities for con_name in ss_ineq_cons: - if "c3" not in con_name and "c5" not in con_name: + should_have_default_priority = ( + all(cname not in con_name for cname in ["c3", "c5", "c12"]) + or "c3_up" in con_name + ) + if should_have_default_priority: self.assertEqual( model_data.separation_priority_order[con_name], - 0, + DEFAULT_SEPARATION_PRIORITY, msg=( f"Separation priority for entry {con_name!r} of second-stage " "inequalities not as expected." @@ -1260,8 +1387,11 @@ def build_simple_test_model_data(self): Build model data object for testing constraint standardization routines. """ - model_data = Bunch() - model_data.config = Bunch() + model_data = ModelData( + original_model=None, + timing=None, + config=Bunch(separation_priority_order=dict()), + ) model_data.working_model = ConcreteModel() model_data.working_model.user_model = m = Block() @@ -1278,6 +1408,7 @@ def build_simple_test_model_data(self): # first-stage equalities m.eq1 = Constraint(expr=m.x1 + log(m.p) == m.q_cert + 1) m.eq2 = Constraint(expr=(1, m.x2, 1)) + m.eq2_unc = Constraint(expr=(m.q, m.x2, m.q)) # second-stage equalities m.eq3 = Constraint(expr=m.x2 * m.q == 1) @@ -1305,6 +1436,7 @@ def build_simple_test_model_data(self): model_data.working_model.original_active_equality_cons = [ m.eq1, m.eq2, + m.eq2_unc, m.eq3, m.eq4, m.eq5, @@ -1326,12 +1458,19 @@ def test_standardize_equality_constraints(self): working_model = model_data.working_model m = working_model.user_model + m.pyros_separation_priority = Suffix() + model_data.config.separation_priority_order[m.eq3.local_name] = 2 + model_data.config.separation_priority_order[m.eq5.local_name] = 3 + m.pyros_separation_priority[m.eq3] = 1 + m.pyros_separation_priority[m.eq4] = 10 + m.pyros_separation_priority[m.eq2_unc] = None + standardize_equality_constraints(model_data) first_stage_eq_cons = working_model.first_stage.dr_independent_equality_cons second_stage_eq_cons = working_model.second_stage.equality_cons - self.assertEqual(len(first_stage_eq_cons), 2) + self.assertEqual(len(first_stage_eq_cons), 3) self.assertEqual(len(second_stage_eq_cons), 4) self.assertFalse(m.eq1.active) @@ -1346,6 +1485,13 @@ def test_standardize_equality_constraints(self): self, new_eq2_con.expr, RangedExpression((1, m.x2, 1), False) ) + self.assertFalse(m.eq2_unc.active) + new_eq2_unc_con = first_stage_eq_cons["eq_con_eq2_unc"] + self.assertTrue(new_eq2_unc_con.active) + assertExpressionsEqual( + self, new_eq2_unc_con.expr, RangedExpression((m.q, m.x2, m.q), False) + ) + self.assertFalse(m.eq3.active) new_eq3_con = second_stage_eq_cons["eq_con_eq3"] self.assertTrue(new_eq3_con.active) @@ -1373,6 +1519,13 @@ def test_standardize_equality_constraints(self): self.assertFalse(m.eq7.active) assertExpressionsEqual(self, m.eq7.expr, m.z2 == 0) + final_priority_dict = model_data.separation_priority_order + self.assertEqual(len(final_priority_dict), 4) + self.assertEqual(final_priority_dict["eq_con_eq3"], 1) + self.assertEqual(final_priority_dict["eq_con_eq4"], 10) + self.assertEqual(final_priority_dict["eq_con_eq5"], 3) + self.assertEqual(final_priority_dict["eq_con_eq6"], DEFAULT_SEPARATION_PRIORITY) + class TestStandardizeActiveObjective(unittest.TestCase): """ @@ -1384,8 +1537,7 @@ def build_simple_test_model_data(self): Build simple model for testing active objective standardization. """ - model_data = Bunch() - model_data.config = Bunch() + model_data = ModelData(original_model=None, timing=None, config=Bunch()) model_data.working_model = ConcreteModel() model_data.working_model.user_model = m = Block() @@ -1429,8 +1581,6 @@ def build_simple_test_model_data(self): model_data.working_model.second_stage = Block() model_data.working_model.second_stage.inequality_cons = Constraint(Any) - model_data.separation_priority_order = dict() - return model_data def test_declare_objective_expressions(self): @@ -1961,8 +2111,14 @@ def setup_test_model_data(self, uncertainty_set=None): Set up simple test model for testing the reformulation routine. """ - model_data = Bunch() - model_data.config = Bunch(uncertainty_set=uncertainty_set or BoxSet([[0, 1]])) + model_data = ModelData( + config=Bunch( + uncertainty_set=uncertainty_set or BoxSet([[0, 1]]), + separation_priority_order=dict(), + ), + original_model=None, + timing=None, + ) model_data.working_model = working_model = ConcreteModel() model_data.working_model.user_model = m = Block() @@ -2000,6 +2156,11 @@ def setup_test_model_data(self, uncertainty_set=None): working_model.second_stage.equality_cons["eq_con_2"] = m.eq_con_2.expr working_model.second_stage.inequality_cons["con"] = m.con.expr + # mock separation priorities added during equality + # constraint standardization + model_data.separation_priority_order["eq_con"] = DEFAULT_SEPARATION_PRIORITY + model_data.separation_priority_order["eq_con_2"] = DEFAULT_SEPARATION_PRIORITY + # deactivate constraints on user model, as these are not # what the reformulation routine actually processes m.eq_con.deactivate() @@ -2082,7 +2243,6 @@ def test_reformulate_nonlinear_state_var_independent_eq_con(self): and recasting of nonlinear constraints to opposing equalities. """ model_data = self.setup_test_model_data() - model_data.separation_priority_order = dict() model_data.config.decision_rule_order = 1 model_data.config.progress_logger = logging.getLogger( @@ -2185,7 +2345,6 @@ def test_reformulate_equality_cons_discrete_set(self): model_data = self.setup_test_model_data( uncertainty_set=DiscreteScenarioSet([[0], [0.7]]) ) - model_data.separation_priority_order = dict() model_data.config.decision_rule_order = 1 model_data.config.progress_logger = logging.getLogger( @@ -2361,6 +2520,10 @@ def build_test_model_data(self): # second-stage equality m.eq4 = Constraint(expr=m.z3 + m.y1 + 5 * m.q2var == m.q) + # duplicate of eq4: we will enforce this only nominally + # through separation priorities + m.eq5 = Constraint(expr=m.z3 + m.y1 + 5 * m.q2var == m.q) + # INEQUALITY CONSTRAINTS # since x1, z1 nonadjustable, LB is first-stage, # but UB second-stage due to uncertain param q @@ -2387,6 +2550,9 @@ def build_test_model_data(self): # to the presence of the uncertain parameter m.ineq6 = Constraint(expr=-m.q2var <= m.x1) + # will be a nominal-only constraint + m.ineq7 = Constraint(expr=m.y3 <= 2 * m.q) + # OBJECTIVE # contains a rich combination of first-stage and second-stage terms m.obj = Objective( @@ -2410,6 +2576,7 @@ def build_test_model_data(self): uncertainty_set=BoxSet([[4, 5], [3, 4], [1, 1]]), uncertain_params=[m.q, m.q2var, m.q_cert], nominal_uncertain_param_vals=[m.q.value, m.q2var.value, m.q_cert.value], + separation_priority_order=dict(), ), ) @@ -2568,9 +2735,15 @@ def test_preprocessor_constraint_partitioning_nonstatic_dr( objective_focus=ObjectiveType[obj_focus], decision_rule_order=dr_order, progress_logger=logger, - separation_priority_order=dict(ineq3=2), + separation_priority_order=dict(ineq3=2, ineq4=3), ) ) + om = model_data.original_model + om.pyros_separation_priority = Suffix() + om.pyros_separation_priority[om.eq5] = None + om.pyros_separation_priority[om.ineq4] = 5 + om.pyros_separation_priority[om.ineq7] = None + preprocess_model_data(model_data, user_var_partitioning) working_model = model_data.working_model @@ -2601,13 +2774,13 @@ def test_preprocessor_constraint_partitioning_nonstatic_dr( self.assertEqual( list(working_model.first_stage.inequality_cons), ( - ["ineq_con_ineq1_lower_bound_con", "ineq_con_ineq2"] + ["ineq_con_ineq1_lower_bound_con", "ineq_con_ineq2", "ineq_con_ineq7"] + (["epigraph_con"] if obj_focus == "nominal" else []) ), ) self.assertEqual( list(working_model.first_stage.dr_independent_equality_cons), - ["eq_con_eq2", "eq_con_eq3"], + ["eq_con_eq2", "eq_con_eq3", "eq_con_eq5"], ) self.assertEqual( list(working_model.first_stage.dr_dependent_equality_cons), @@ -2733,6 +2906,9 @@ def test_preprocessor_constraint_partitioning_nonstatic_dr( ss.inequality_cons["ineq_con_ineq6_lower_bound_con"].expr, -m.x1 <= -(-1 * working_model.temp_uncertain_params[1]), ) + assertExpressionsEqual( + self, fs.inequality_cons["ineq_con_ineq7"].expr, m.y3 <= 2 * m.q + ) assertExpressionsEqual( self, fs.dr_independent_equality_cons["eq_con_eq2"].expr, m.x1 - m.z1 == 0 @@ -2745,21 +2921,35 @@ def test_preprocessor_constraint_partitioning_nonstatic_dr( if dr_order < 2: # due to coefficient matching, this should have been deleted self.assertNotIn("eq_con_eq1", ss.equality_cons) + assertExpressionsEqual( + self, + fs.dr_independent_equality_cons["eq_con_eq5"].expr, + m.z3 + m.y1 + 5 * working_model.temp_uncertain_params[1] == m.q, + ) # user model block should have no active constraints self.assertFalse(list(m.component_data_objects(Constraint, active=True))) - # check separation priorities + # check separation priorities are as expected + self.assertEqual( + list(model_data.separation_priority_order.keys()), + list(ss.inequality_cons.keys()), + ) + final_priority_dict = model_data.separation_priority_order + self.assertEqual(final_priority_dict["ineq_con_ineq3_lower_bound_con"], 2) + self.assertEqual(final_priority_dict["ineq_con_ineq3_upper_bound_con"], 2) + self.assertEqual(final_priority_dict["ineq_con_ineq4_lower_bound_con"], 5) for con_name, order in model_data.separation_priority_order.items(): - expected_order = 2 if "ineq3" in con_name else 0 - self.assertEqual( - order, - expected_order, - msg=( - "Separation priority order for second-stage inequality " - f"{con_name!r} not as expected." - ), - ) + if "ineq3" not in con_name and "ineq4" not in con_name: + self.assertEqual( + order, + DEFAULT_SEPARATION_PRIORITY, + msg=( + "Separation priority order for second-stage inequality " + f"{con_name!r} not as expected." + ), + ) + self.assertFalse(ublk.pyros_separation_priority.active) @parameterized.expand( [["static", 0, True], ["affine", 1, False], ["quadratic", 2, False]] @@ -2777,7 +2967,7 @@ def test_preprocessor_coefficient_matching( objective_focus=ObjectiveType.worst_case, decision_rule_order=dr_order, progress_logger=logger, - separation_priority_order=dict(), + separation_priority_order=dict(eq1=1), ) ) @@ -2785,15 +2975,23 @@ def test_preprocessor_coefficient_matching( # due to the coefficient matching constraints derived # from bounds on z5 robust_infeasible = preprocess_model_data(model_data, user_var_partitioning) - self.assertIsInstance(robust_infeasible, bool) - self.assertEqual(robust_infeasible, expected_robust_infeas) - # check the coefficient matching constraint expressions working_model = model_data.working_model m = model_data.working_model.user_model fs = working_model.first_stage dr_dependent_fs_eqs = working_model.first_stage.dr_dependent_equality_cons ss_ineqs = working_model.second_stage.inequality_cons + + self.assertIsInstance(robust_infeasible, bool) + self.assertEqual(robust_infeasible, expected_robust_infeas) + if not expected_robust_infeas: + # all equality constraints were processed, + # so only inequality constraint names should appear in the + # priority dict + self.assertEqual( + list(model_data.separation_priority_order.keys()), list(ss_ineqs.keys()) + ) + if config.decision_rule_order == 1: # check the constraint expressions of eq1 and z5 bound assertExpressionsEqual( @@ -2847,6 +3045,11 @@ def test_preprocessor_coefficient_matching( m.q * (m.z3 + m.x2 * m.q_cert) <= 0.0, ) + # check separation priority properly accounted for + sep_priority_dict = model_data.separation_priority_order + self.assertEqual(sep_priority_dict["reform_upper_bound_from_eq_con_eq1"], 1) + self.assertEqual(sep_priority_dict["reform_lower_bound_from_eq_con_eq1"], 1) + # check coefficient matching constraint expressions assertExpressionsEqual( self, @@ -2942,6 +3145,40 @@ def test_preprocessor_objective_standardization(self, name, dr_order): ), ) + @parameterized.expand([["nominal"], ["worst_case"]]) + def test_preprocessor_sep_priorities_suffix_finder_none(self, obj_focus): + """ + Test preprocessor resolves separation priorities as expected + when an active separation priority Suffix component + has a custom value mapped to `None`. + """ + model_data, user_var_partitioning = self.build_test_model_data() + config = model_data.config + config.update( + dict( + objective_focus=ObjectiveType[obj_focus], + decision_rule_order=1, + progress_logger=logger, + separation_priority_order=dict(), + ) + ) + model_data.original_model.pyros_separation_priority = Suffix() + model_data.original_model.pyros_separation_priority[None] = 10 + preprocess_model_data(model_data, user_var_partitioning) + ss_ineq_cons = model_data.working_model.second_stage.inequality_cons + self.assertEqual( + list(model_data.separation_priority_order.keys()), list(ss_ineq_cons.keys()) + ) + for con_idx in ss_ineq_cons.keys(): + if con_idx != "epigraph_con": + self.assertEqual(model_data.separation_priority_order[con_idx], 10) + else: + # custom prioritization of epigraph constraint is ignored + self.assertEqual( + model_data.separation_priority_order[con_idx], + DEFAULT_SEPARATION_PRIORITY, + ) + @parameterized.expand([["nominal"], ["worst_case"]]) def test_preprocessor_log_model_statistics_affine_dr(self, obj_focus): """ @@ -2955,7 +3192,7 @@ def test_preprocessor_log_model_statistics_affine_dr(self, obj_focus): objective_focus=ObjectiveType[obj_focus], decision_rule_order=1, progress_logger=logger, - separation_priority_order=dict(), + separation_priority_order=dict(eq5=None, ineq7=None), ) ) preprocess_model_data(model_data, user_var_partitioning) @@ -2971,14 +3208,14 @@ def test_preprocessor_log_model_statistics_affine_dr(self, obj_focus): State variables : 2 (1 adj.) Decision rule variables : 6 Number of uncertain parameters : 3 (2 eff.) - Number of constraints : 26 - Equality constraints : 11 + Number of constraints : 28 + Equality constraints : 12 Coefficient matching constraints : 6 - Other first-stage equations : 2 + Other first-stage equations : 3 Second-stage equations : 1 Decision rule equations : 2 - Inequality constraints : 15 - First-stage inequalities : {3 if obj_focus == 'nominal' else 2} + Inequality constraints : 16 + First-stage inequalities : {4 if obj_focus == 'nominal' else 3} Second-stage inequalities : {12 if obj_focus == 'nominal' else 13} """ ) @@ -3007,7 +3244,7 @@ def test_preprocessor_log_model_statistics_quadratic_dr(self, obj_focus): objective_focus=ObjectiveType[obj_focus], decision_rule_order=2, progress_logger=logger, - separation_priority_order=dict(), + separation_priority_order=dict(eq5=None, ineq7=None), ) ) preprocess_model_data(model_data, user_var_partitioning) @@ -3023,14 +3260,14 @@ def test_preprocessor_log_model_statistics_quadratic_dr(self, obj_focus): State variables : 2 (1 adj.) Decision rule variables : 12 Number of uncertain parameters : 3 (2 eff.) - Number of constraints : 28 - Equality constraints : 11 + Number of constraints : 30 + Equality constraints : 12 Coefficient matching constraints : 6 - Other first-stage equations : 2 + Other first-stage equations : 3 Second-stage equations : 1 Decision rule equations : 2 - Inequality constraints : 17 - First-stage inequalities : {3 if obj_focus == 'nominal' else 2} + Inequality constraints : 18 + First-stage inequalities : {4 if obj_focus == 'nominal' else 3} Second-stage inequalities : {14 if obj_focus == 'nominal' else 15} """ ) diff --git a/pyomo/contrib/pyros/tests/test_separation.py b/pyomo/contrib/pyros/tests/test_separation.py index 14fcb6b9478..548551726c4 100644 --- a/pyomo/contrib/pyros/tests/test_separation.py +++ b/pyomo/contrib/pyros/tests/test_separation.py @@ -305,6 +305,7 @@ def test_construct_separation_problem_uncertain_factor_param_components(self): class TestGroupSecondStageIneqConsByPriority(unittest.TestCase): def test_group_ss_ineq_constraints_by_priority(self): model_data = build_simple_model_data() + model_data.config.progress_logger = logging.getLogger(__name__) separation_model = construct_separation_problem(model_data) # build mock separation data-like object @@ -312,6 +313,7 @@ def test_group_ss_ineq_constraints_by_priority(self): separation_data = Bunch( separation_model=separation_model, separation_priority_order=model_data.separation_priority_order, + config=model_data.config, ) priority_groups = group_ss_ineq_constraints_by_priority(separation_data) diff --git a/pyomo/contrib/pyros/util.py b/pyomo/contrib/pyros/util.py index 257d2529bf5..3fee93ecc48 100644 --- a/pyomo/contrib/pyros/util.py +++ b/pyomo/contrib/pyros/util.py @@ -25,6 +25,7 @@ from pyomo.common.collections import ComponentMap, ComponentSet from pyomo.common.dependencies import scipy as sp +from pyomo.common.flags import NOTSET from pyomo.common.errors import ApplicationError, InvalidValueError from pyomo.common.log import Preformatted from pyomo.common.modeling import unique_component_name @@ -46,6 +47,7 @@ VarData, value, ) +from pyomo.core.base.suffix import SuffixFinder from pyomo.core.expr.numeric_expr import SumExpression from pyomo.core.expr.numvalue import native_types from pyomo.core.expr.visitor import ( @@ -75,6 +77,7 @@ TIC_TOC_SOLVE_TIME_ATTR = "pyros_tic_toc_time" DEFAULT_LOGGER_NAME = "pyomo.contrib.pyros" DEFAULT_SEPARATION_PRIORITY = 0 +BYPASSING_SEPARATION_PRIORITY = None class TimingData: @@ -1083,6 +1086,11 @@ class ModelData: separation_priority_order : dict Mapping from constraint names to separation priority values. + separation_priority_suffix_finder : SuffixFinder + Object for resolving active Suffix components + added to the model by the user as a means of + prioritizing separation problems mapped to + second-stage inequality constraints. """ def __init__(self, original_model, config, timing): @@ -1092,6 +1100,9 @@ def __init__(self, original_model, config, timing): self.separation_priority_order = dict() # working model will be addressed by preprocessing self.working_model = None + self.separation_priority_suffix_finder = SuffixFinder( + name="pyros_separation_priority", default=NOTSET + ) def preprocess(self, user_var_partitioning): """ @@ -1106,6 +1117,46 @@ def preprocess(self, user_var_partitioning): """ return preprocess_model_data(self, user_var_partitioning) + def get_user_separation_priority(self, component_data, component_data_name): + """ + Infer user specification for the separation priority/priorities + of the second-stage inequality constraint/constraints derived + from a given component data attribute of the working model. + + Parameters + ---------- + component_data : ComponentData + Component data from which the inequality constraints + are meant to be derived. + component_data_name : str + Name of the component data object as it is expected + to appear in ``self.config.separation_priority_order``. + + Returns + ------- + numeric type or None + Priority of the derived constraint(s). + + Notes + ----- + The separation priorities for the constraints derived from + ``component_data`` are inferred from either the + active ``Suffix`` components of ``self.working_model`` + with name 'pyros_separation_priority' + or from ``self.config.separation_priority``. + Priorities specified through the active ``Suffix`` components + take precedence over priorities specified through + ``self.config.separation_priority_order``. + Moreover, priorities are inferred from ``Suffix`` components + using the Pyomo ``SuffixFinder``. + """ + priority = self.separation_priority_suffix_finder.find(component_data) + if priority is NOTSET: + priority = self.config.separation_priority_order.get( + component_data_name, DEFAULT_SEPARATION_PRIORITY + ) + return priority + def setup_quadratic_expression_visitor( wrt, subexpression_cache=None, var_map=None, var_order=None, sorter=None @@ -1626,6 +1677,10 @@ def turn_nonadjustable_var_bounds_to_constraints(model_data): Consequently, all constraints added to the working model in this method are considered second-stage constraints. + Uncertain bounds that have been assigned a separation priority + of None are not reformulated, and are subsequently enforced + subject to only the nominal uncertain parameter realization. + Parameters ---------- model_data : model data object @@ -1640,11 +1695,18 @@ def turn_nonadjustable_var_bounds_to_constraints(model_data): var_name = var.getname( relative_to=working_model.user_model, fully_qualified=True ) + var_bound_sep_priority = model_data.get_user_separation_priority( + component_data=var, component_data_name=var_name + ) for btype, bound in declared_bound_triple._asdict().items(): is_bound_uncertain = bound is not None and ( ComponentSet(identify_mutable_parameters(bound)) & uncertain_params_set ) - if is_bound_uncertain: + is_bound_second_stage = ( + is_bound_uncertain + and var_bound_sep_priority is not BYPASSING_SEPARATION_PRIORITY + ) + if is_bound_second_stage: new_con_expr = create_bound_constraint_expr(var, bound, btype) new_con_name = f"var_{var_name}_uncertain_{btype}_bound_con" remove_var_declared_bound(var, btype) @@ -1656,15 +1718,9 @@ def turn_nonadjustable_var_bounds_to_constraints(model_data): working_model.second_stage.inequality_cons[new_con_name] = ( new_con_expr ) - # can't specify custom priorities for variable bounds - model_data.separation_priority_order[new_con_name] = ( - DEFAULT_SEPARATION_PRIORITY - ) - - # for subsequent developments: return a mapping - # from each variable to the corresponding binding constraints? - # we will add this as needed when changes are made to - # the interface for separation priority ordering + model_data.separation_priority_order[new_con_name] = ( + var_bound_sep_priority + ) def turn_adjustable_var_bounds_to_constraints(model_data): @@ -1680,6 +1736,11 @@ def turn_adjustable_var_bounds_to_constraints(model_data): Since these constraints depend on adjustable variables, they are taken to be (effective) second-stage constraints. + Bounds that have been assigned a separation priority + of None are reformulated to first-stage constraints, + and are subsequently enforced subject to only the + nominal uncertain parameter realization. + Parameters ---------- model_data : model data object @@ -1702,32 +1763,35 @@ def turn_adjustable_var_bounds_to_constraints(model_data): ("certain", cert_bound_triple), ("uncertain", uncert_bound_triple), ) + var_bound_sep_priority = model_data.get_user_separation_priority( + component_data=var, component_data_name=var_name + ) + is_bound_second_stage = ( + var_bound_sep_priority is not BYPASSING_SEPARATION_PRIORITY + ) + if is_bound_second_stage: + ineq_con_component = working_model.second_stage.inequality_cons + eq_con_component = working_model.second_stage.equality_cons + else: + ineq_con_component = working_model.first_stage.inequality_cons + eq_con_component = working_model.first_stage.dr_independent_equality_cons for certainty_desc, bound_triple in cert_uncert_bound_zip: for btype, bound in bound_triple._asdict().items(): if bound is not None: new_con_name = f"var_{var_name}_{certainty_desc}_{btype}_bound_con" new_con_expr = create_bound_constraint_expr(var, bound, btype) if btype == BoundType.EQ: - working_model.second_stage.equality_cons[new_con_name] = ( - new_con_expr - ) + eq_con_component[new_con_name] = new_con_expr else: - working_model.second_stage.inequality_cons[new_con_name] = ( - new_con_expr - ) - # no custom separation priorities for Var - # bound constraints + ineq_con_component[new_con_name] = new_con_expr + + if is_bound_second_stage: model_data.separation_priority_order[new_con_name] = ( - DEFAULT_SEPARATION_PRIORITY + var_bound_sep_priority ) remove_all_var_bounds(var) - # for subsequent developments: return a mapping - # from each variable to the corresponding binding constraints? - # we will add this as needed when changes are made to - # the interface for separation priority ordering - def _replace_vars_in_component_exprs(block, substitution_map, ctype): """ @@ -1913,7 +1977,6 @@ def standardize_inequality_constraints(model_data): model_data : model data object Main model data object, containing the working model. """ - config = model_data.config working_model = model_data.working_model uncertain_params_set = ComponentSet(working_model.effective_uncertain_params) adjustable_vars_set = ComponentSet( @@ -1930,8 +1993,13 @@ def standardize_inequality_constraints(model_data): con_rel_name = con.getname( relative_to=working_model.user_model, fully_qualified=True ) - - if uncertain_params_in_con_expr | adjustable_vars_in_con_body: + con_sep_priority = model_data.get_user_separation_priority( + component_data=con, component_data_name=con_rel_name + ) + is_con_potentially_second_stage = ( + uncertain_params_in_con_expr | adjustable_vars_in_con_body + ) and con_sep_priority is not BYPASSING_SEPARATION_PRIORITY + if is_con_potentially_second_stage: con_bounds_triple = rearrange_bound_pair_to_triple( lower_bound=con.lower, upper_bound=con.upper ) @@ -1969,9 +2037,7 @@ def standardize_inequality_constraints(model_data): ) # account for user-specified priority specifications model_data.separation_priority_order[new_con_name] = ( - config.separation_priority_order.get( - con_rel_name, DEFAULT_SEPARATION_PRIORITY - ) + con_sep_priority ) else: # we do not want to modify the arrangement of @@ -2024,14 +2090,20 @@ def standardize_equality_constraints(model_data): con_rel_name = con.getname( relative_to=working_model.user_model, fully_qualified=True ) - if uncertain_params_in_con_expr | adjustable_vars_in_con_body: - working_model.second_stage.equality_cons[f"eq_con_{con_rel_name}"] = ( + con_sep_priority = model_data.get_user_separation_priority( + component_data=con, component_data_name=con_rel_name + ) + new_con_name = f"eq_con_{con_rel_name}" + is_con_second_stage = ( + uncertain_params_in_con_expr | adjustable_vars_in_con_body + ) and con_sep_priority is not BYPASSING_SEPARATION_PRIORITY + if is_con_second_stage: + working_model.second_stage.equality_cons[new_con_name] = con.expr + model_data.separation_priority_order[new_con_name] = con_sep_priority + else: + working_model.first_stage.dr_independent_equality_cons[new_con_name] = ( con.expr ) - else: - working_model.first_stage.dr_independent_equality_cons[ - f"eq_con_{con_rel_name}" - ] = con.expr # definitely don't want active duplicate con.deactivate() @@ -2529,7 +2601,7 @@ def _reformulate_eq_con_continuous_uncertainty( working_model.second_stage.inequality_cons[new_con_name] = std_con_expr # no custom priorities specified model_data.separation_priority_order[new_con_name] = ( - DEFAULT_SEPARATION_PRIORITY + model_data.separation_priority_order[ss_eq_con_index] ) else: polynomial_repn_coeffs = ( @@ -2724,6 +2796,7 @@ def reformulate_state_var_independent_eq_cons(model_data): ) if robust_infeasible: break + del model_data.separation_priority_order[con_idx] # we no longer need these auxiliary components working_model.del_component(temp_param_vars) @@ -2830,6 +2903,10 @@ def preprocess_model_data(model_data, user_var_partitioning): ) robust_infeasible = reformulate_state_var_independent_eq_cons(model_data) + # we are done looking for separation priorities + for priority_sfx in model_data.separation_priority_suffix_finder.all_suffixes: + priority_sfx.deactivate() + return robust_infeasible