1# ___________________________________________________________________________ 2# 3# Pyomo: Python Optimization Modeling Objects 4# Copyright 2017 National Technology and Engineering Solutions of Sandia, LLC 5# Under the terms of Contract DE-NA0003525 with National Technology and 6# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain 7# rights in this software. 8# This software is distributed under the 3-clause BSD License. 9# ___________________________________________________________________________ 10 11import pyomo.core.expr.current as EXPR 12from pyomo.gdp import GDP_Error, Disjunction 13from pyomo.gdp.disjunct import _DisjunctData, Disjunct 14 15from pyomo.core.base.component import _ComponentBase 16 17from pyomo.core import Block, TraversalStrategy 18from pyomo.opt import TerminationCondition, SolverStatus 19from weakref import ref as weakref_ref 20import logging 21 22logger = logging.getLogger('pyomo.gdp') 23 24_acceptable_termination_conditions = set([ 25 TerminationCondition.optimal, 26 TerminationCondition.globallyOptimal, 27 TerminationCondition.locallyOptimal, 28]) 29_infeasible_termination_conditions = set([ 30 TerminationCondition.infeasible, 31 TerminationCondition.invalidProblem, 32]) 33 34 35class NORMAL(object): pass 36class INFEASIBLE(object): pass 37class NONOPTIMAL(object): pass 38 39def verify_successful_solve(results): 40 status = results.solver.status 41 term = results.solver.termination_condition 42 43 if status == SolverStatus.ok and term in _acceptable_termination_conditions: 44 return NORMAL 45 elif term in _infeasible_termination_conditions: 46 return INFEASIBLE 47 else: 48 return NONOPTIMAL 49 50 51def clone_without_expression_components(expr, substitute=None): 52 """A function that is used to clone an expression. 53 54 Cloning is roughly equivalent to calling ``copy.deepcopy``. 55 However, the :attr:`clone_leaves` argument can be used to 56 clone only interior (i.e. non-leaf) nodes in the expression 57 tree. Note that named expression objects are treated as 58 leaves when :attr:`clone_leaves` is :const:`True`, and hence 59 those subexpressions are not cloned. 60 61 This function uses a non-recursive 62 logic, which makes it more scalable than the logic in 63 ``copy.deepcopy``. 64 65 Args: 66 expr: The expression that will be cloned. 67 substitute (dict): A dictionary mapping object ids to 68 objects. This dictionary has the same semantics as 69 the memo object used with ``copy.deepcopy``. Defaults 70 to None, which indicates that no user-defined 71 dictionary is used. 72 73 Returns: 74 The cloned expression. 75 """ 76 if substitute is None: 77 substitute = {} 78 # 79 visitor = EXPR.ExpressionReplacementVisitor(substitute=substitute, 80 remove_named_expressions=True) 81 return visitor.dfs_postorder_stack(expr) 82 83def preprocess_targets(targets): 84 preprocessed_targets = [] 85 for t in targets: 86 if t.ctype is Disjunction: 87 if t.is_indexed(): 88 for disjunction in t.values(): 89 for disj in disjunction.disjuncts: 90 preprocessed_targets.append(disj) 91 else: 92 for disj in t.disjuncts: 93 preprocessed_targets.append(disj) 94 # now we are safe to put the disjunction, and if the target was 95 # anything else, then we don't need to worry because disjuncts 96 # are declared before disjunctions they appear in 97 preprocessed_targets.append(t) 98 return preprocessed_targets 99 100def target_list(x): 101 if isinstance(x, _ComponentBase): 102 return [ x ] 103 elif hasattr(x, '__iter__'): 104 ans = [] 105 for i in x: 106 if isinstance(i, _ComponentBase): 107 ans.append(i) 108 else: 109 raise ValueError( 110 "Expected Component or list of Components." 111 "\n\tReceived %s" % (type(i),)) 112 return ans 113 else: 114 raise ValueError( 115 "Expected Component or list of Components." 116 "\n\tReceived %s" % (type(x),)) 117 118# [ESJ 07/09/2019 Should this be a more general utility function elsewhere? I'm 119# putting it here for now so that all the gdp transformations can use it. 120# Returns True if child is a node or leaf in the tree rooted at parent, False 121# otherwise. Accepts list of known components in the tree and updates this list 122# to enhance performance in future calls. Note that both child and parent must 123# be blocks! 124def is_child_of(parent, child, knownBlocks=None): 125 # Note: we can get away with a dictionary and not ComponentMap because we 126 # will only store Blocks (or their ilk), and Blocks are hashable (only 127 # derivatives of NumericValue are not hashable) 128 if knownBlocks is None: 129 knownBlocks = {} 130 tmp = set() 131 node = child 132 while True: 133 known = knownBlocks.get(node) 134 if known: 135 knownBlocks.update({c: True for c in tmp}) 136 return True 137 if known is not None and not known: 138 knownBlocks.update({c: False for c in tmp}) 139 return False 140 if node is parent: 141 knownBlocks.update({c: True for c in tmp}) 142 return True 143 if node is None: 144 knownBlocks.update({c: False for c in tmp}) 145 return False 146 147 tmp.add(node) 148 container = node.parent_component() 149 if container is node: 150 node = node.parent_block() 151 else: 152 node = container 153 154def get_src_disjunction(xor_constraint): 155 """Return the Disjunction corresponding to xor_constraint 156 157 Parameters 158 ---------- 159 xor_constraint: Constraint, which must be the logical constraint 160 (located on the transformation block) of some 161 Disjunction 162 """ 163 # NOTE: This is indeed a linear search through the Disjunctions on the 164 # model. I am leaving it this way on the assumption that asking XOR 165 # constraints for their Disjunction is not going to be a common 166 # question. If we ever need efficiency then we should store a reverse 167 # map from the XOR constraint to the Disjunction on the transformation 168 # block while we do the transformation. And then this method could query 169 # that map. 170 m = xor_constraint.model() 171 for disjunction in m.component_data_objects(Disjunction, 172 descend_into=(Block, Disjunct)): 173 if disjunction._algebraic_constraint: 174 if disjunction._algebraic_constraint() is xor_constraint: 175 return disjunction 176 raise GDP_Error("It appears that '%s' is not an XOR or OR constraint " 177 "resulting from transforming a Disjunction." 178 % xor_constraint.name) 179 180def get_src_disjunct(transBlock): 181 """Return the Disjunct object whose transformed components are on 182 transBlock. 183 184 Parameters 185 ---------- 186 transBlock: _BlockData which is in the relaxedDisjuncts IndexedBlock 187 on a transformation block. 188 """ 189 if not hasattr(transBlock, "_srcDisjunct") or \ 190 type(transBlock._srcDisjunct) is not weakref_ref: 191 raise GDP_Error("Block '%s' doesn't appear to be a transformation " 192 "block for a disjunct. No source disjunct found." 193 % transBlock.name) 194 return transBlock._srcDisjunct() 195 196def get_src_constraint(transformedConstraint): 197 """Return the original Constraint whose transformed counterpart is 198 transformedConstraint 199 200 Parameters 201 ---------- 202 transformedConstraint: Constraint, which must be a component on one of 203 the BlockDatas in the relaxedDisjuncts Block of 204 a transformation block 205 """ 206 transBlock = transformedConstraint.parent_block() 207 # This should be our block, so if it's not, the user messed up and gave 208 # us the wrong thing. If they happen to also have a _constraintMap then 209 # the world is really against us. 210 if not hasattr(transBlock, "_constraintMap"): 211 raise GDP_Error("Constraint '%s' is not a transformed constraint" 212 % transformedConstraint.name) 213 # if something goes wrong here, it's a bug in the mappings. 214 return transBlock._constraintMap['srcConstraints'][transformedConstraint] 215 216def _find_parent_disjunct(constraint): 217 # traverse up until we find the disjunct this constraint lives on 218 parent_disjunct = constraint.parent_block() 219 while not isinstance(parent_disjunct, _DisjunctData): 220 if parent_disjunct is None: 221 raise GDP_Error( 222 "Constraint '%s' is not on a disjunct and so was not " 223 "transformed" % constraint.name) 224 parent_disjunct = parent_disjunct.parent_block() 225 226 return parent_disjunct 227 228def _get_constraint_transBlock(constraint): 229 parent_disjunct = _find_parent_disjunct(constraint) 230 # we know from _find_parent_disjunct that parent_disjunct is a Disjunct, 231 # so the below is OK 232 transBlock = parent_disjunct._transformation_block 233 if transBlock is None: 234 raise GDP_Error("Constraint '%s' is on a disjunct which has not been " 235 "transformed" % constraint.name) 236 # if it's not None, it's the weakref we wanted. 237 transBlock = transBlock() 238 239 return transBlock 240 241def get_transformed_constraints(srcConstraint): 242 """Return the transformed version of srcConstraint 243 244 Parameters 245 ---------- 246 srcConstraint: ScalarConstraint or _ConstraintData, which must be in 247 the subtree of a transformed Disjunct 248 """ 249 if srcConstraint.is_indexed(): 250 raise GDP_Error("Argument to get_transformed_constraint should be " 251 "a ScalarConstraint or _ConstraintData. (If you " 252 "want the container for all transformed constraints " 253 "from an IndexedDisjunction, this is the parent " 254 "component of a transformed constraint originating " 255 "from any of its _ComponentDatas.)") 256 transBlock = _get_constraint_transBlock(srcConstraint) 257 try: 258 return transBlock._constraintMap['transformedConstraints'][srcConstraint] 259 except: 260 logger.error("Constraint '%s' has not been transformed." 261 % srcConstraint.name) 262 raise 263 264def _warn_for_active_disjunction(disjunction, disjunct, NAME_BUFFER): 265 # this should only have gotten called if the disjunction is active 266 assert disjunction.active 267 problemdisj = disjunction 268 if disjunction.is_indexed(): 269 for i in sorted(disjunction.keys()): 270 if disjunction[i].active: 271 # a _DisjunctionData is active, we will yell about 272 # it specifically. 273 problemdisj = disjunction[i] 274 break 275 276 parentblock = problemdisj.parent_block() 277 # the disjunction should only have been active if it wasn't transformed 278 assert problemdisj.algebraic_constraint is None 279 _probDisjName = problemdisj.getname( 280 fully_qualified=True, name_buffer=NAME_BUFFER) 281 _disjName = disjunct.getname(fully_qualified=True, name_buffer=NAME_BUFFER) 282 raise GDP_Error("Found untransformed disjunction '%s' in disjunct '%s'! " 283 "The disjunction must be transformed before the " 284 "disjunct. If you are using targets, put the " 285 "disjunction before the disjunct in the list." 286 % (_probDisjName, _disjName)) 287 288def _warn_for_active_disjunct(innerdisjunct, outerdisjunct, NAME_BUFFER): 289 assert innerdisjunct.active 290 problemdisj = innerdisjunct 291 if innerdisjunct.is_indexed(): 292 for i in sorted(innerdisjunct.keys()): 293 if innerdisjunct[i].active: 294 # This shouldn't be true, we will complain about it. 295 problemdisj = innerdisjunct[i] 296 break 297 298 raise GDP_Error("Found active disjunct '{0}' in disjunct '{1}'! Either {0} " 299 "is not in a disjunction or the disjunction it is in " 300 "has not been transformed. {0} needs to be deactivated " 301 "or its disjunction transformed before {1} can be " 302 "transformed.".format( 303 problemdisj.getname( 304 fully_qualified=True, name_buffer = NAME_BUFFER), 305 outerdisjunct.getname( 306 fully_qualified=True, 307 name_buffer=NAME_BUFFER))) 308 309 310def _warn_for_active_logical_constraint(logical_statement, disjunct, NAME_BUFFER): 311 # this should only have gotten called if the logical constraint is active 312 assert logical_statement.active 313 problem_statement = logical_statement 314 if logical_statement.is_indexed(): 315 for i in logical_statement: 316 if logical_statement[i].active: 317 # a _LogicalConstraintData is active, we will yell about 318 # it specifically. 319 problem_statement = logical_statement[i] 320 break 321 # None of the _LogicalConstraintDatas were actually active. We 322 # are OK and we can deactivate the container. 323 else: 324 logical_statement.deactivate() 325 return 326 # the logical constraint should only have been active if it wasn't transformed 327 _probStatementName = problem_statement.getname( 328 fully_qualified=True, name_buffer=NAME_BUFFER) 329 raise GDP_Error("Found untransformed logical constraint %s in disjunct %s! " 330 "The logical constraint must be transformed before the " 331 "disjunct. Use the logical_to_linear transformation." 332 % (_probStatementName, disjunct.name)) 333 334 335def check_model_algebraic(instance): 336 """Checks if there are any active Disjuncts or Disjunctions reachable via 337 active Blocks. If there are not, it returns True. If there are, it issues 338 a warning detailing where in the model there are remaining non-algebraic 339 components, and returns False. 340 341 Parameters 342 ---------- 343 instance: a Model or Block 344 """ 345 disjunction_set = {i for i in instance.component_data_objects( 346 Disjunction, descend_into=(Block, Disjunct), active=None)} 347 active_disjunction_set = {i for i in instance.component_data_objects( 348 Disjunction, descend_into=(Block, Disjunct), active=True)} 349 disjuncts_in_disjunctions = set() 350 for i in disjunction_set: 351 disjuncts_in_disjunctions.update(i.disjuncts) 352 disjuncts_in_active_disjunctions = set() 353 for i in active_disjunction_set: 354 disjuncts_in_active_disjunctions.update(i.disjuncts) 355 356 for disjunct in instance.component_data_objects( 357 Disjunct, descend_into=(Block,), 358 descent_order=TraversalStrategy.PostfixDFS): 359 # check if it's relaxed 360 if disjunct.transformation_block is not None: 361 continue 362 # It's not transformed, check if we should complain 363 elif disjunct.active and _disjunct_not_fixed_true(disjunct) and \ 364 _disjunct_on_active_block(disjunct): 365 # If someone thinks they've transformed the whole instance, but 366 # there is still an active Disjunct on the model, we will warn 367 # them. In the future this should be the writers' job.) 368 if disjunct not in disjuncts_in_disjunctions: 369 logger.warning('Disjunct "%s" is currently active, ' 370 'but was not found in any Disjunctions. ' 371 'This is generally an error as the model ' 372 'has not been fully relaxed to a ' 373 'pure algebraic form.' % (disjunct.name,)) 374 return False 375 elif disjunct not in disjuncts_in_active_disjunctions: 376 logger.warning('Disjunct "%s" is currently active. While ' 377 'it participates in a Disjunction, ' 378 'that Disjunction is currently deactivated. ' 379 'This is generally an error as the ' 380 'model has not been fully relaxed to a pure ' 381 'algebraic form. Did you deactivate ' 382 'the Disjunction without addressing the ' 383 'individual Disjuncts?' % (disjunct.name,)) 384 return False 385 else: 386 logger.warning('Disjunct "%s" is currently active. It must be ' 387 'transformed or deactivated before solving the ' 388 'model.' % (disjunct.name,)) 389 return False 390 # We didn't find anything bad. 391 return True 392 393def _disjunct_not_fixed_true(disjunct): 394 # Return true if the disjunct indicator variable is not fixed to True 395 return not (disjunct.indicator_var.fixed and 396 disjunct.indicator_var.value) 397 398def _disjunct_on_active_block(disjunct): 399 # Check first to make sure that the disjunct is not a descendent of an 400 # inactive Block or fixed and deactivated Disjunct, before raising a 401 # warning. 402 parent_block = disjunct.parent_block() 403 while parent_block is not None: 404 # deactivated Block 405 if parent_block.ctype is Block and not parent_block.active: 406 return False 407 # properly deactivated Disjunct 408 elif (parent_block.ctype is Disjunct and not parent_block.active 409 and parent_block.indicator_var.value == False 410 and parent_block.indicator_var.fixed): 411 return False 412 else: 413 # Step up one level in the hierarchy 414 parent_block = parent_block.parent_block() 415 continue 416 return True 417