-
Notifications
You must be signed in to change notification settings - Fork 577
KKT Transform Code Refactor #3881
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Changes from 43 commits
c9c1283
5797b93
32de8cd
9ec1e58
a4537f0
57c94bb
6cf46c6
fce4632
a41c5ee
a9149b1
d62f32c
74638a7
115019d
cc4ebde
4f3b7f7
cff252f
c10d5a7
f5ea2b0
7cf29ff
cfa0555
e961546
9df1401
e04ac57
2de227e
08c82f2
10c7e10
7583fd4
c0fe0d0
77e926a
8eec4a2
42f9b49
30614f1
dc16852
4a994c0
d7f91d9
deb48d7
253dede
09e08eb
ceb0fda
e4c8a56
009fa86
915979f
55370d3
4457093
b5117d0
e8aca52
1ce03dc
1bec2cc
3c0dab0
ca735e5
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -21,4 +21,5 @@ | |
| scaling, | ||
| logical_to_linear, | ||
| lp_dual, | ||
| kkt, | ||
| ) | ||
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,251 @@ | ||
| # ____________________________________________________________________________________ | ||
| # | ||
| # Pyomo: Python Optimization Modeling Objects | ||
| # Copyright (c) 2008-2026 National Technology and Engineering Solutions of Sandia, LLC | ||
| # Under the terms of Contract DE-NA0003525 with National Technology and Engineering | ||
| # Solutions of Sandia, LLC, the U.S. Government retains certain rights in this | ||
| # software. This software is distributed under the 3-clause BSD License. | ||
| # ____________________________________________________________________________________ | ||
|
|
||
|
|
||
| from pyomo.common.autoslots import AutoSlots | ||
| from pyomo.common.collections import ComponentMap, ComponentSet | ||
| from pyomo.common.config import ConfigDict, ConfigValue | ||
| from pyomo.core import ( | ||
| Block, | ||
| Constraint, | ||
| ConstraintList, | ||
| Expression, | ||
| NonNegativeReals, | ||
| Objective, | ||
| TransformationFactory, | ||
| Var, | ||
| VarList, | ||
| ) | ||
| from pyomo.core.expr.calculus.diff_with_pyomo import reverse_sd | ||
| from pyomo.core.expr.visitor import identify_variables | ||
| from pyomo.mpec import ComplementarityList, complements | ||
| from pyomo.util.config_domains import ComponentDataSet | ||
|
|
||
|
|
||
| class _KKTReformulationData(AutoSlots.Mixin): | ||
| __slots__ = ("obj_dual_map", "dual_obj_map") | ||
|
|
||
| def __init__(self): | ||
| self.obj_dual_map = ComponentMap() | ||
| self.dual_obj_map = ComponentMap() | ||
|
|
||
|
|
||
| Block.register_private_data_initializer(_KKTReformulationData) | ||
|
|
||
|
|
||
| @TransformationFactory.register( | ||
| 'core.kkt', 'Generate KKT reformulation of the given model' | ||
| ) | ||
| class NonLinearProgrammingKKT: | ||
| CONFIG = ConfigDict("core.kkt") | ||
| CONFIG.declare( | ||
| 'kkt_block_name', | ||
| ConfigValue( | ||
| default='kkt', | ||
| doc=""" | ||
| Name of the block on which the kkt variables and constraints will be stored. | ||
| """, | ||
| ), | ||
| ) | ||
| CONFIG.declare( | ||
| 'parametrize_wrt', | ||
| ConfigValue( | ||
| default=[], | ||
| domain=ComponentDataSet(Var), | ||
| description='Vars to treat as data for the purposes of generating KKT reformulation', | ||
| doc=""" | ||
| Optional list of Vars to be treated as data while generating the KKT reformulation. | ||
| """, | ||
| ), | ||
| ) | ||
|
|
||
| def apply_to(self, model, **kwds): | ||
| """ | ||
| Reformulate model with KKT conditions. | ||
| """ | ||
| config = self.CONFIG(kwds.pop('options', {})) | ||
| config.set_value(kwds) | ||
|
|
||
| if hasattr(model, config.kkt_block_name): | ||
| raise ValueError( | ||
| "model already has an attribute with the " | ||
| f"specified kkt_block_name: '{config.kkt_block_name}'" | ||
| ) | ||
|
|
||
| # we should check that all vars the user fixed are included | ||
| # in parametrize_wrt | ||
|
ChrisLaliwala marked this conversation as resolved.
Outdated
|
||
| params = config.parametrize_wrt | ||
|
|
||
| kkt_block = Block(concrete=True) | ||
| kkt_block.parametrize_wrt = params | ||
| self._reformulate(model, kkt_block, params) | ||
| model.add_component(config.kkt_block_name, kkt_block) | ||
| return model | ||
|
|
||
| def _reformulate(self, model, kkt_block, params): | ||
| # initialize | ||
| info = model.private_data() | ||
| lagrangean = 0 | ||
| all_vars_set = ComponentSet() | ||
|
|
||
| # collect the active Objectives | ||
| active_objs = list( | ||
| model.component_data_objects(Objective, active=True, descend_into=True) | ||
| ) | ||
| if len(active_objs) != 1: | ||
| raise ValueError( | ||
| f"model must have only one active objective; found {len(active_objs)}" | ||
|
ChrisLaliwala marked this conversation as resolved.
Outdated
|
||
| ) | ||
| # collect vars from active objective | ||
| obj = active_objs[0] | ||
| all_vars_set.update(identify_variables(obj.expr, include_fixed=True)) | ||
| lagrangean += obj.sense * obj.expr | ||
|
|
||
| # list of equality multipliers | ||
| kkt_block.gamma = VarList() | ||
| # list of inequality multipliers | ||
| kkt_block.alpha = VarList(domain=NonNegativeReals) | ||
| # define inequality complements | ||
| kkt_block.complements = ComplementarityList() | ||
|
|
||
| for con in model.component_data_objects( | ||
| Constraint, descend_into=True, active=True | ||
| ): | ||
| lower, body, upper = con.to_bounded_expression() | ||
|
|
||
| # collect variables in constraint | ||
| for expr in (lower, body, upper): | ||
| if expr is None: | ||
| continue | ||
| all_vars_set.update(identify_variables(expr=expr, include_fixed=True)) | ||
|
|
||
| if con.equality: | ||
| gamma_i = kkt_block.gamma.add() | ||
| lagrangean += (upper - body) * gamma_i | ||
| info.obj_dual_map[con] = gamma_i | ||
| info.dual_obj_map[gamma_i] = con | ||
|
|
||
| else: | ||
| alpha_l = None | ||
| if lower is not None: | ||
| alpha_l = kkt_block.alpha.add() | ||
| con_expr = lower - body | ||
| lagrangean += con_expr * alpha_l | ||
| kkt_block.complements.add(complements(alpha_l >= 0, con_expr <= 0)) | ||
| info.dual_obj_map[alpha_l] = con | ||
|
|
||
| alpha_u = None | ||
| if upper is not None: | ||
| alpha_u = kkt_block.alpha.add() | ||
| con_expr = body - upper | ||
| lagrangean += con_expr * alpha_u | ||
| kkt_block.complements.add(complements(alpha_u >= 0, con_expr <= 0)) | ||
| info.dual_obj_map[alpha_u] = con | ||
|
|
||
| info.obj_dual_map[con] = (alpha_l, alpha_u) | ||
|
|
||
| fixed_vars = ComponentSet(v for v in all_vars_set if v.is_fixed()) | ||
| var_set = ComponentSet(all_vars_set) | ||
| var_set -= fixed_vars | ||
|
|
||
| # do error checking on parametrize_wrt | ||
| missing = fixed_vars - params | ||
| if missing: | ||
| raise ValueError("All fixed variables must be included in parametrize_wrt.") | ||
|
emma58 marked this conversation as resolved.
Outdated
|
||
|
|
||
| if not params <= all_vars_set: | ||
| raise ValueError( | ||
| "A variable passed in parametrize_wrt does not exist on an " | ||
|
ChrisLaliwala marked this conversation as resolved.
Outdated
|
||
| "active constraint or objective within the model." | ||
| ) | ||
|
|
||
| var_set = var_set - params | ||
| for var in var_set: | ||
| alpha_l = None | ||
| if var.has_lb(): | ||
| alpha_l = kkt_block.alpha.add() | ||
| con_expr = var.lb - var | ||
| lagrangean += con_expr * alpha_l | ||
| kkt_block.complements.add(complements(alpha_l >= 0, con_expr <= 0)) | ||
| info.dual_obj_map[alpha_l] = var | ||
|
|
||
| alpha_u = None | ||
| if var.has_ub(): | ||
| alpha_u = kkt_block.alpha.add() | ||
| con_expr = var - var.ub | ||
| lagrangean += con_expr * alpha_u | ||
| kkt_block.complements.add(complements(alpha_u >= 0, con_expr <= 0)) | ||
| info.dual_obj_map[alpha_u] = var | ||
|
|
||
| info.obj_dual_map[var] = (alpha_l, alpha_u) | ||
|
|
||
| kkt_block.lagrangean = Expression(expr=lagrangean) | ||
|
|
||
| # enforce stationarity conditions | ||
| deriv_lagrangean = reverse_sd(kkt_block.lagrangean.expr) | ||
| kkt_block.stationarity_conditions = ConstraintList() | ||
| for var in var_set: | ||
| kkt_block.stationarity_conditions.add(deriv_lagrangean[var] == 0) | ||
|
|
||
| active_objs[0].deactivate() | ||
|
|
||
| def get_object_from_multiplier(self, model, multiplier_var): | ||
|
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I would prefer this be named
Member
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I have to admit that I would prefer
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I was afraid you would say that! :P OK, let's leave it then. |
||
| """ | ||
| Return the constraint corresponding to a KKT multiplier variable. If the | ||
| multiplier corresponds to an inequality formed by a variable bound, the variable | ||
| is returned. | ||
|
|
||
| Parameters | ||
| ---------- | ||
| model: ConcreteModel | ||
| The model on which the kkt transformation was applied | ||
| multiplier_var: Var | ||
| A KKT multiplier created by the transformation. | ||
|
|
||
| Returns | ||
| ------- | ||
| Object | ||
| - Constraint object | ||
| - Variable | ||
| """ | ||
|
|
||
| info = model.private_data() | ||
| if multiplier_var in info.dual_obj_map: | ||
| return info.dual_obj_map[multiplier_var] | ||
| raise ValueError( | ||
| f"The KKT multiplier: {multiplier_var.name}, does not exist on {model.name}." | ||
| ) | ||
|
|
||
| def get_multiplier_from_object(self, model, component): | ||
|
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Same naming thought as above |
||
| """ | ||
| Return the multiplier for the object. If the object is a normal constraint, a single | ||
| multiplier is returned. If the object is a ranged constraint or a variable, a tuple | ||
| containing the lower and upper bound multipliers is returned. | ||
|
|
||
| Parameters | ||
| ---------- | ||
| model: ConcreteModel | ||
| The model on which the kkt transformation was applied to | ||
|
ChrisLaliwala marked this conversation as resolved.
Outdated
|
||
| component: Constraint or Variable | ||
|
|
||
| Returns | ||
| ------- | ||
| VarData | tuple[VarData | None, VarData | None] | ||
| The KKT multiplier(s) corresponding to the component. | ||
| For ranged constraints/variables, returns (lb_mult, ub_mult), | ||
| where an entry is 'None' if that bound doesn't exist. | ||
| """ | ||
|
|
||
| info = model.private_data() | ||
| if component in info.obj_dual_map: | ||
| return info.obj_dual_map[component] | ||
| raise ValueError( | ||
| f"The component '{component.name}' either does not exist on " | ||
| f"'{model.name}', or is not associated with a multiplier." | ||
| ) | ||
Uh oh!
There was an error while loading. Please reload this page.