[PEP] Solving Blocks with variables that do not live on the Block

[michaelbynum]

I would like to be able to solve Pyomo Blocks with variables that do not live on the Block. I will describe my use case, but there are others.

I have a structured MINLP. In fact, it looks like a binary tree of Blocks. I have a transformation to relax this MINLP. The the root block has a constraint with a bilinear term (for example), x*y, then I will replace that bilinear term with an auxiliary variable, aux, and add constraints for the McCormick envelopes relating aux to x and y. If this same bilinear term, x*y appears in other constraints, I want to replace it with the same auxiliary variable aux. I do not want to create a new variable for the same term. Now, suppose I encounter the same bilinear term in a constraint on a child Block. I already have the auxiliary variable on the root Block. It is not easy to move the auxiliary variable so that it has a new parent. I also have a need to solve individual child Blocks.

Another example I am aware of arises in GDP transformations. Maybe @jsiirola can elaborate.

[carldlaird]

I see the use case you describe. Two comments:

1) Not that I think this is better, but would it make sense to use a transformation to flatten the model first? Then this kind of issue would go away? (I suspect this would bring a host of other issues, not the least of which is the required name mangling when flattening.)

2) Is this (and other use cases like this) implying that variables are not really tied to their parent?

[emma58]

I think that besides the naming issues, flattening the model also obscures the structure which was likely why you were using Blocks in the first place. I don't think there is any reason to assume that you are going to stop messing with a model after just one solve. If you have to flatten to do that solve, you've just made it much harder to make further changes. And I think if individual Blocks on the original model were solvable before the transformation, it seems reasonable to expect them to be after it as well. Granted, transformations could be designed so that this was still true, but in the GDP transformations at least, we are using Blocks so that the transformed model contains a completely intact version of the original model (possibly with some components deactivated). Everything the transformation adds is on a new Block it creates. This helps map the original model to what the transformation did and makes the names of components created by transformations predictable.

[qtothec]

I'll add that flattening also obscures the structure that the solvers I am involved with rely upon.

[michaelbynum]

The only reason I can think of to require that variables live on the block bening solved is for model validation (which certainly depends on how you conceptually define a model). However, model validation is not really what is currently happening. Consider the following:

>>> b = pe.Block(concrete=True)
>>> b.x = pe.Var()
>>> b.b = pe.Block(concrete=True)
>>> b.b.obj = pe.Objective(expr=b.x**2)
>>> b.b.c = pe.Constraint(expr=b.x >= 1)
>>> opt = pe.SolverFactory('ipopt')
>>> res = opt.solve(b)
>>> res = opt.solve(b.b)
Traceback (most recent call last):
...
KeyError: 90652546328

Here, b can be solved, but b.b cannot. However, in both cases, constraints are defined which use variables that do not live on the same block.

[michaelbynum]

@carldlaird

1. Flattening the model would solve this problem. However, the block structure was actually induced on purpose and is necessary.

2. I’m not sure. There might be two questions here. First, how is an optimization model defined? Second, what are the software design issues? For design reasons, it may be necessary to “tie” a variable to it’s parent - I’m not sure.

[michaelbynum]

Here is another use case. Suppose I want to model a stochastic programming problem, but I want to use the same Var objects (for the first stage variables) in all of the scenarios rather than using nonanticipativity constraints. I can’t think of a good reason for this right now, but it is certainly a valid way to model stochastic programs. I could put the first stage variables on the root block, and use a block for each of the scenarios. I could then use the first stage variables in each of the scenario sub-blocks, but I would not be able to solve the scenario sub-blocks independently.

[jwatsonnm]

Regarding the above case, an alternative would be to: (1) define a block containing all first stage variables, (2) reference that block in both the scenario blocks and the master model block, (3) define the master optimization objective by referencing cost terms in the individual scenario blocks. Then you at least have encapsulation. But you have other issues with respect to non-tree block reference structures. I bring up this alternative because I really don't think violating encapsulation is a good idea (e.g., are we going to allow variable names of the type "child_block."..".first_stage_var?). If violating encapsulation is a good idea, I assert that we should first ensure that supporting a block DAG in model specification isn't a cleaner or at least equally viable alternative.

[carldlaird]

I like the way you are leaning here JP. Basically, if violating encapsulation is a good idea, then this is an indicator that we do not have the correct structure for grouping constraints (e.g., you said a DAG would be better than a tree).

[michaelbynum]

My primary question (which I don't know the answer to) is whether or not variables should be involved in encapsulation. Why do we add variables to models/blocks? Is it simply for convenient access so that they may be used in other components (e.g., constraints, objectives, disjuncts, etc.), or is it actually for model definition?

If we want to enforce encapsulation, then should we throw an error when a constraint uses a variable that lives elsewhere in the model? We do not currently do this (as illustrated by the example above).

A block DAG is an interesting idea. This would require blocks to have more than one parent though, right?

[ghackebeil]

I like this suggestion as well.

I would argue, however, that there is a simple way to handle Michael's SP example without any changes to the code. Assuming the master variables are on the top-level model, if you want to solve the scenarios individually, simply do the following:

# deactivate all scenarios
for scenario in model.scenarios:
    scenario.deactivate()

# solve each scenario
for scenario in model.scenarios:
    scenario.activate()
    solve(model)
    scenario.deactivate()

[carldlaird]

@ghackebeil Can you be more specific? Which suggestion do you like?

[michaelbynum]

@ghackebeil I agree that your solution works for the two stage problem. It would be more complex in the multi-stage case (which is more similar to my original use case). Additionally, what if I don't want to include the constraints on the root block? I'm thinking of a bender's type algorithm...

[ghackebeil]

That was meant for suggestion given by @jwatsonnm.

@michaelbynum: then you could deactivate the constraints on the root block.

[michaelbynum]

I am not trying to argue that there are not alternative ways to handle this issue. However, I am trying to ask a more fundamental question. Should this be possible? If the answer is no, I am okay with that. I think it is worth discussing though.

[ghackebeil]

I agree that it is worth discussing. I like the idea about allowing objects to be placed on more than one block / in more than one container. It would allow you to put together a temporary model from a collections of previously declared components. Perhaps these special containers wouldn't change the parent pointer on the children that they store. Component labels could then be generated in two different ways:

1. Bottom up, resulting in possible collisions

2. Top down, no collisions

   This reminds me of how temporary models are put together in AMPL.

[jwatsonnm]

I personally believe variables - and actually all components - should be involved in the encapsulation. So in my opinion it's a model definition issue. If we don't support an all-or-nothing approach regarding components and encapsulation, then I feel we can't claim an object-oriented mathematical programming model. Which is a defining feature of Pyomo. In general, I assert every component should "live" somewhere in a clean design. Or, perhaps more correctly, requiring such will lead / guide users toward clean designs.

Correct in that the DAG model would have multiple parents in general for blocks. If you follow encapsulation strictly, then this would not matter too much, i.e., you wouldn't reference "up" if you were encapsulating.

I would be a proponent of optionally throwing exceptions - in a model checker capability sense - if encapsulation is violated.

[carldlaird]

Another option that "fixes" this conundrum is if we supported variable aliasing - not the same as a reference variable, but true aliasing. The ability to say, I have a variable here and another variable here, but they are actually the same physical variable. Then you would write your blocks "encapsulated", and if you wanted one variable at the root, you would use aliases instead of equality constraints.

Not sure this captures all use cases, but it would capture many.

[jwatsonnm]

Following up on most recent comment by @ghackebeil above: A nice use case that we can't support right now in Pyomo are, for lack of a better term, common blocks. For example, global data that can be referenced by multiple blocks in a hierarchical model.

[michaelbynum]

I agree that it is worth discussing. I like the idea about allowing objects to be placed on more than one block / in more than one container. It would allow you to put together a temporary model from a collections of previously declared components. Perhaps these special containers wouldn't change the parent pointer on the children that they store. Component labels could then be generated in two different ways:

1. Bottom up, resulting in possible collisions
2. Top down, no collisions

This reminds me of how temporary models are put together in AMPL.

I really like this idea.

[carldlaird]

@jwatsonnm Not sure exactly what you mean here. If you mean "data", then I think this is completely supported today (if you violate encapsulation and reference data on another block).

[jwatsonnm]

@carldlaird I'd like to not violate encapsulation and reference data on a parent block in the hierarchy.

[carldlaird]

@michaelbynum @ghackebeil Regarding the idea of components on more than one block - this would require some thought - I think that these could be two separate data structures. There is one (with encapsulation) for the purposes of ownership, but another data structure for describing what to include in a particular model that you want the solver to solve (a different container of objectives, constraints, etc.)

[carldlaird]

@jwatsonnm Can you be more specific? I don't understand what you mean by not violating encapsulation and also referencing data on a parent block in the hierarchy?

[carldlaird]

@jwatsonnm I have advocated for some time that we include an additional argument to rules to pass additional data when the component is created. This would work for blocks or other components. Maybe this is what you mean? Something like:

def con_rule(m, i, extra_passed_data):
      return m.x[I] == extra_passed_data.my_desired_x[i]
my_block.con = Constraint(set_I, rule=con_rule, extra_data=jps_global_data)

It looks like encapsulation is violated in the final model structure that is created, but the "component definition" is reusable since this is explicitly passed in. Components could be included in that extra data object as well.

def con_rule(m, i, extra_passed_data):
      return m.x[I] == extra_passed_data['global_x'][i]
my_block.con = Constraint(set_I, rule=con_rule, extra_data={'global_x': parent.x})

instead of something like:

def con_rule(m, i):
      return m.x[I] == m.parent().x[i]
my_block.con = Constraint(set_I, rule=con_rule)

[ghackebeil]

I’m not sure I understand why that needs to be accessed through an argument to the rule. It’s already available in global scope, correct?

On May 28, 2019, at 12:30 PM, carldlaird notifications@github.com wrote:

@jwatsonnm I have advocated for some time that we include an additional argument to rules to pass additional data when the component is created. This would work for blocks or other components. Maybe this is what you mean? Something like:

def con_rule(m, i, extra_passed_data): return m.x[I] == extra_passed_data.my_desired_x[i] my_block.con = Constraint(set_I, rule=con_rule, extra_data=jps_global_data) — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[carldlaird]

@ghackebeil Doesn't this depend on when and where the rule is fired?

[ghackebeil]

Never mind. It depends on where the rules are defined.

On May 28, 2019, at 12:33 PM, carldlaird notifications@github.com wrote:

It looks like encapsulation is violated in the final model structure that is created, but the "component definition" is reusable since this is explicitly passed in. Components could be included in that extra data object as well.

def con_rule(m, i, extra_passed_data): return m.x[I] == extra_passed_data['global_x'][i] my_block.con = Constraint(set_I, rule=con_rule, extra_data={'global_x': parent.x}) instead of something like:

def con_rule(m, i): return m.x[I] == m.parent().x[i] my_block.con = Constraint(set_I, rule=con_rule) — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[jsiirola]

I am in favor of allowing Var objects to live outside the tree of active blocks that defines the current call to solve(). I am even potentially amenable to allowing Vars that do not live in the same model tree (either on a separate tree or on no tree at all) ... I think it is bad practice to do it, but allowing it probably follows from my arguments below.

First, GDP currently has a need for this, as the indicator variable that is assigned to the Disjunct is a part of the relaxed model, however, the Disjunct is deactivated after the constraints defined in it have all been relaxed (onto another block). Moving the Var is generally undesirable, as the user will not have an easy way to "find" it again after the transformation. Our current solution is to go through the model and reactivate Disjuncts after they all have been relaxed, plus deactivating all the embedded original constraints. This process is a hack at best, and has proven to be very fragile, especially for nested disjunctions (which we do see in the wild). That said, I will admit that part of the problem is that Pyomo does not differentiate between Boolean and Binary variables. If we did, then we would need to "relax" the Boolean var into a Binary var as part of the transformation, and then pull the resulting state back to the Boolean var after the solve. However, implementing that will involve standing up a new logical expression system, which will be an involved PEP in itself.

PAO has some interesting use cases for variable declaration as well, where the inner problem is allowed to reference variable (values, in this case) from the outer problem.

Now, my argument for why we should treat Variables differently from Constraints arises from my proposal of what it means to be "a model". I propose the following definition:

A model is the set of Active component objects reachable through Active objects from a single reference Block.

Most components in Pyomo are ActiveComponent derivatives. This includes Blocks, Constraints, etc. It is notable the list of objects that are not ActiveComponents in Pyomo: Var, Param, Set (and Set-like things), and Suffix. Why are they not derived from ActiveComponent? Mostly, it is because it is not obvious what it means to "deactivate" them. For example, what is a deactivated Set? The argument for Var is a little trickier: per my definition above (which I admit is open to debate), a deactivated Var would not be "part of the model". This is different from "fixing" the variable, where the value is still part of the model. The closest thing that I can come up with would be like "removing the column" (in the sense of branch-and-price) from the A-matrix for MIPs. While this makes sense for MIPs, it is less clear to me what would/should happen for nonlinear constraints. In my mind, this ambiguity is the biggest reason we don't allow deactivating variables. OK - so if that makes sense for a single Var, what about a Var on a deactivated sub-block? That Block is not part of the model (again, per my definition). But, if a Constraint in the active portion of the problem (say, the top block) refers to the Var, what should the solver do? I think omitting the Var is probably not correct (for the same arguments as deactivating a Var), so in my mind, the path of "least surprise" is to pass it to the solver. Finally, if you are going to allow a Var on a deactivated block, then it is a small step to allow references to a Var on another Block outside the current active subtree passed to solve().

Second argument: with our current situation (where Variables are implicitly required to exist in the solve()'s current active subtree), we have the "same" information - the variables in the model - stored twice: once through the Vars reachable on the current active subtree, and once by calling identify_variables() on the active constraints in the current active subtree. In general, this violates "good database design", where every piece of information is stored in exactly one way. For database designers, this rule arose because storing a piece of information twice was a form of caching, and ensuring the validity of the data store meant that the consumer had to make sure that the two locations were always kept in sync (i.e., solving the cache invalidation problem). We have exactly the same problem, as the two lists of variables (through iterating blocks and through iterating constraints) are already expected and allowed to disagree: we allow variables to exist in the subtree but not in any constraints. We do not change their value when the model is solved, either: we simply pretend they didn't exist (except for marking them "stale"). Enforcing "correct" models in the current situation (where Variables must be in the subtree) is actually a somewhat expensive operation, as we have to build (and store) both lists and check that they agree. This seems to me to be "extra" work that we don't need to do (especially when we are trying to improve the overall performance of the solution process).

Finally, I think we should draw careful distinction between containment and encapsulation (in the OO sense). Blocks are containers that hold modeling components and can help to organize them into useful representations. However, I do not think that the tree structure that Blocks enable is the only useful structure for organizing a model. While you can use a tree structure to represent a set of independent objects, I think we will continue to find modeling needs that aren't conveniently assembled by composing sub-trees. I hit this early-on in the Grid work (SCUC with transmission switching), and IDAES sees this now. Indeed, IDAES is heavily using References, which can be used (abused) to have some modeling components to appear in multiple places in the tree. (This complicates things and makes it very easy to emit redundant constraints). While I think we should strongly encourage our users to leverage encapsulation concepts when they build models, I do not think Pyomo should require it.

I really like @jwatsonnm's idea of moving the "encapsulation test" to a model checker-like functionality that a user can run occasionally during development.

I am interested in supporting alternative representations (like DAG structures), but I suspect that it will be complicated (I am not entirely sure how to define encapsulation at that point), and will take a longer design discussion with some interesting motivating examples.

[michaelbynum]

I would like to revisit this issue.

It is clear that it is very important to do some error checking when handing a model to a solver. However, not allowing solution of blocks with variables that live elsewhere is also clearly problematic (especially for complex workflows). I would like to start documenting options along with their advantages and disadvantages. Please feel free to help me flush these out. If I missed one of the suggestions above, I apologize.

1. Provide functions and other tools which do model validation. Users would be required to call these functions explicitly if they want model validation. Calls to solve would, by default, use any variables found in the ActiveComponents encountered. Calls to solve would have options to enforce model validation if desired.

   • Advantages
     • Very flexible
     • We could still issue warnings by default if constraints are found which use variables not defined on the block handed to the solver. However, this might get annoying in some cases.
   • Disadvantages
     • Users who do not utilize the model validation tools may miss bugs in their code
     • If we adopted this approach, many examples would need to be updated to utilize the model validation tools so that users know they exist and to get users into the habit of using them.

2. Add support for aliases. If a user wants to solve a block with variables defined elsewhere, they must put an alias to the variable on the block.

   • Advantages
     • Default model validation
   • Disadvantages
     • This may impact Pyomo iterators such as component_data_objects. For example, we probably do not want the same variable to appear twice in list(b.component_data_objects(Var)).

3. Do not allow the solution of blocks with variables defined elsewhere.

   • Advantages
     • default model validation
   • Disadvantages
     • not flexible

4. If a user wants to solve a block with variables defined elsewhere, then the "extra" variables must be handed to the solver with a special keyword argument.

   • Advantages
     • default model validation
   • Disadvantages
     • Although this is more flexible than option (3), it is still not flexible, especially for complex workflows with transformations, third-party packages, etc.

5. Allow the solution of blocks as long as the variables used in the constraints on the block live somewhere on the root model

6. Add flags to blocks to indicate whether or not they “know” all of their variables.

[michaelbynum]

Note that options (1) and (4) are essentially the same with different defaults.

[qtothec]

There's a looser option to (4), which is "trust me, I know what I'm doing"

[jsiirola]

One thing that I noticed when working on the LP writer rewrite is that "default model validation" in options 2-4 does actually require additional work and may be a bit of a misnomer; it may be better to think of it as "enforced model validation".

I am also beginning to think that there is a scale between "anything goes" for variable use and the strict "all variables must be on active subblocks". For example, I am beginning to think that Variables that are not assigned to any block may in fact indicate a modeling error.

If we are collecting proposals, I will toss in a variant of @michaelbynum's (5):

Writers will collect all variables used in active constraints reachable from the block being solved through a tree of active blocks (the "active tree") and verify that they all share a common model(). In addition, any fixed variables appearing in the constraints or attached to the active tree will have their bounds verified; infeasible fixed variables will be logged as errors. Any free variables in the active tree not appearing in constraints will either be sent to the solver as free variables and/or regularized (set to the valid value closest to 0) during solution storage.

[qtothec]

I just realized that automatically suppressing warnings is also not too difficult, even without passing a flag all the way to where the warning would be raised: https://github.com/Pyomo/pyomo/blob/d78034e7123a8ccf6b3801a9f53c6707fb3f0778/pyomo/contrib/gdpopt/util.py#L45

Relevant for certain options.

[michaelbynum]

@qtothec That is very relevant and useful. Thanks for sharing.

[michaelbynum]

I am mostly adding this comment for future reference (instead of relying on my memory). Consider the relax_integrality transformation. This transformation takes a block or model and iterates over the variables using component_data_objects. This implicitly assumes that all of the variables used in the block or model can be found through component_data_objects. While this is not enforced in any way (the transformation will happily run even if component_data_objects does not yield anything), it is assumed. I believe a great deal of Pyomo-based code relies on this undocumented implication. If we decide that models are defined by the set of active components on them, then I think we need to deprecate iterating over variables using component_data_objects.