Running a Concrete model multiple times in the same script

[Christopher Hagmann]

I have a script that in essence looks like

def create_model(data):
   model = ConcreteModel()
   
   #initialize Param with data

def solve_model(data):

    model = create_model(data)

    #solve model in script

if __name__ == '__main__':

    

    instance1 = solve_model(data1)

    instance2 = solve_model(data2)

I could attached the actual script but it is really long. Solving the first instance works great but the second one always raises an error. I can swap the two lines and the first will still solve and the second not so I know it isn't a issue with the data input. The error is posted below.

Traceback (most recent call last):

  File "BigMModel.py", line 725, in <module>

    instance, T = solve_big_m_model(gap=.02)

  File "/home/cdhagmann/Documents/GitHub/WITP/Function_Module.py", line 252, in _inner

    Results = func(*args, **kwargs)

  File "BigMModel.py", line 467, in solve_big_m_model

    model = big_m_model(PUTAWAY, PICKING)

  File "BigMModel.py", line 250, in big_m_model

    model.constraint1 = Constraint(expr=summation(model.theta_put) == 1)

  File "/opt/coopr/src/coopr.pyomo/coopr/pyomo/base/numvalue.py", line 448, in __eq__

    return generate_relational_expression(_eq, self, other)

  File "/opt/coopr/src/coopr.pyomo/coopr/pyomo/base/expr.py", line 1411, in generate_relational_expression

    raise TypeError(chainedInequalityErrorMessage())

TypeError: Nonconstant relational expression used in an unexpected Boolean context.

The inequality expression:

    0.0  <=  MHE_Cost[i2]

contains non-constant terms (variables) appearing in a Boolean context, e.g.:

    if expression <= 5:

This is generally invalid.  If you want to obtain the Boolean value of

the expression based on the current variable values, explicitly evaluate

the expression, e.g.:

    if value(expression) <= 5:

or

    if value(expression <= 5):

It seems like it is trying to reuse the first model but I don't know why or how. How should I be structuring this?

[Christopher Hagmann]

So in trying things out to see if I could debug this on my own I moved `model.constraint1` after the the second constraint and now it works. Why would it being added to the model first or second matter in solving the model?

[Watson, Jean-Paul]

There is no reason that I can conceive of for why that would be the case. Can you post the complete model, so we can investigate?

jpw

--
You received this message because you are subscribed to the Google Groups "Pyomo Forum" group.
To unsubscribe from this group and stop receiving emails from it, send an email to pyomo-forum...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Christopher Hagmann]

I am now getting a different bug that is having the same effect where I cannot seem to be able to run the same model twice consistantly. I am running Pyomo 4.0, but I was having the same problem on Pyomo 3.5 which is why I had upgraded in the first place. Code is at https://github.com/cdhagmann/BigM

I have a concrete model that I cannot run twice on "large" models. It will run twice on small data sets and it will run once on larger data sets, but not twice. Below is the output is you were to run python BigMModel.py:

['i1', 'i2', 'i3', 'i4', 'i5', 'i6'] ['j1', 'j2', 'j3', 'j4', 'j5', 'j6']
Selected Tech10: $42,316.35 [17.86%] (19.11 seconds)
19.11 seconds
['i1', 'i2', 'i3', 'i4', 'i5', 'i6'] ['j1', 'j2', 'j3', 'j4', 'j5', 'j6']
ERROR: Rule failed when generating expression for constraint BigM_MHE_LOWER with index i1:
AttributeError: 'NoneType' object has no attribute 'cname'
ERROR: Constructing component 'BigM_MHE_LOWER' from data=None failed:
AttributeError: 'NoneType' object has no attribute 'cname'
Traceback (most recent call last):
File "BigMModel.py", line 721, in <module>
instance, T = solve_big_m_model(gap=.2)
File "BigMModel.py", line 461, in solve_big_m_model
instance = big_m_model(PUTAWAY, PICKING)
File "BigMModel.py", line 251, in big_m_model
model.BigM_MHE_LOWER = Constraint(model.PUTAWAY, rule=BigM_MHE_LOWER)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/block.py", line 396, in __setattr__
self.add_component(name, val)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/block.py", line 676, in add_component
val.construct(data)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/constraint.py", line 340, in construct
tmp = apply_indexed_rule(self, _self_rule, _self_parent, val)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/misc.py", line 59, in apply_indexed_rule
return rule(model, index)
File "BigMModel.py", line 249, in BigM_MHE_LOWER
return -model.M_MHE * (1 - model.theta_put[i]) <= expr
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/numvalue.py", line 434, in __le__
return generate_relational_expression(_le, self, other)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/expr.py", line 1411, in generate_relational_expression
raise TypeError(chainedInequalityErrorMessage())
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/expr.py", line 66, in chainedInequalityErrorMessage
generate_relational_expression.chainedInequality.to_string(buf)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/expr.py", line 491, in to_string
precedence=_my_precedence )
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/expr.py", line 771, in to_string
precedence=_sub_precedence )
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/component.py", line 488, in to_string
ostream.write(self.__str__())
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/component.py", line 482, in __str__
return self.cname(True)
File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/component.py", line 499, in cname
base = c.cname(fully_qualified, name_buffer)
AttributeError: 'NoneType' object has no attribute 'cname'

[Gabriel Hackebeil]

I was able to run your script without issue using Pyomo trunk, so maybe the problem has been resolved.

Gabe

[Christopher Hagmann]

Gabe,

I just installed the trunk this morning through svn, and I am consistently getting this error on the instance that I posted. Could this have anything to do with this being on a virtual machine?

[Gabe Hackebeil]

I suspect this is a reference counting issue specific to the version of Python you are running, but we're both on Python 2.7.x so that's not entirely certain. Now it looks like you're bumping into an exception while trying generate the "real" error message, so it's hard to say what the real problem is. Can you add a print statement for what you return from that constraint rule? Also send whatever version information you can about your Python.

Gabe

[Siirola, John D]

There are two errors.  The first is the chained inequality logic throwing an exception.  See the stack trace:

 

File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/expr.py", line 1411, in generate_relational_expression

raise TypeError(chainedInequalityErrorMessage())

 

This is usually do to models that perform relational (<= or >= ) comparisons with model components in if clauses without enclosing the relevant component in a call to “value()”.

 

Unfortunately, when generating the exception message, we hit a situation where one of the expressions is invalid.  This leads to another exception while trying to generate the string name for the component data object in the offending relational expression: the component data object has no parent component! 

 

File "/usr/local/lib/python2.7/dist-packages/Pyomo-4.0-py2.7.egg/pyomo/core/base/component.py", line 499, in cname

base = c.cname(fully_qualified, name_buffer)

AttributeError: 'NoneType' object has no attribute 'cname'

 

SO, here’s my guess:

1)      Chris is performing a comparison somewhere in his logic between model construction and reconstructing the model that includes model components without calling value().  This triggers some of the inequality chaining logic (this is deep in the bowels of the expression generation system).

2)      To pull off inequality chaining, we hold on to a reference to the first inequality that we see (so we know how to properly construct the chained inequality).

3)      When Chris goes to reconstruct the model, the original component is removed, and garbage collection collects the component, but because the chained inequality system holds a reference to the expression (which in turn refers to the component data), the component data is not collected.  Now we have a component data object with no parent component object.

4)      When we get around to constructing the next relational expression, the inequality chaining system gets upset and attempts to raise an exception.  This then leads to the stack trace Chris sent.

 

Alas, I am not entirely sure what the best way to trap this and return a meaningful error message.

 

Now, as to why Chris sees this and Gabe does not:  My guess is that the machines have different amounts of RAM, and so the two interpreters feel different pressure as to how frequently they need to garbage collect: if garbage collection doesn’t fire, then the second exception won’t be raised, and Chris would see a more meaningful error from the chained inequality generation system. 

 

At least that’s my guess … and there is certainly a significant probability that I’m wrong.

 

john

[Gabriel Hackebeil]

Chris,

I do see number of lines in your model suggesting that John’s analysis is correct. E.g.,

if model.tb <= t <= model.ty:

   …

and

[… if model.tb <= t <= model.ty]

Since model.tb and model.ty are Params, it is likely that you are (a) not getting the result you would expect (always True) and (b) causing the chained inequality error somewhere after.

Gabe

[Christopher Hagmann]

I am only doing if statements with model parameters. Do they need some sort of value function? I do agree it has to do something with the garbage collection, because I just found out that If I store the first model in memory the second model works just fine.

...

[Siirola, John D]

You need value() with all mutable Params (if you think about it long enough, fixed Vars and mutable Params are basically the same thing).  I think you also need value() with singleton params:  it is possible that you don’t (and I have always wanted to relax that restriction), but I don’t remember implementing that logic..

 

It is technically safe to omit the value() for indexed immutable Params, but Gabe will argue that you should always use value() for safety and consistency.  Even though it is (marginally) slower, he is probably right.  ;)

--

[Gabe Hackebeil]

Since you are working with concrete data, it's advisable to drop the use of Param altogether (except in cases where you need something mutable). This would be the fastest alternative and likely require little to no changes to the constraints in your model. If you do stick with Params, John is correct in that I do advise wrapping them in value() anytime you use them outside of a constraint expression.

Gabe 

[Christopher Hagmann]

I changed all of the singleton parameter (particularly the ones used in `if` logic) to just plain integers and I am still getting the same error.

[Markus Hartikainen]

I am having a similar problem in running the problem

from pyomo.environ import *
import itertools
from pyomo.opt import SolverFactory
import numpy as np
polytope1=list({0,1,2})
polytope2=list({3,4,5})
points = np.random.rand(6,3)
epsilon = 1e-3

model = ConcreteModel()

model.polytope1_coefficient = Var(polytope1,bounds=(0.0,1.0),domain=NonNegativeReals)
model.polytope2_coefficient = Var(polytope2,bounds=(0.0,1.0),domain=NonNegativeReals)
model.polytope1 = Constraint(expr=summation(model.polytope1_coefficient) == 1, doc="Coefficients of polytope 1 must sum up to one")
model.polytope2 = Constraint(expr=summation(model.polytope2_coefficient) == 1, doc="Coefficients of polytope 2 must sum up to one")
model.objective = Objective(expr=max([sum([model.polytope1_coefficient[corner]*points[corner][objective] for corner in polytope1])- sum([model.polytope2_coefficient[corner]*points[corner][objective] for corner in polytope2]) for objective in range(len(points[0]))]), sense=minimize, doc='Define objective function')
opt = SolverFactory("glpk")
opt.solve(model)
print (model.objective.expr() < -epsilon)

Doesn not seem to run twice in a row. I get an error message

TypeError                                 Traceback (most recent call last)
<ipython-input-2-55f91f9db70c> in <module>()
     12 model.polytope1_coefficient = Var(polytope1,bounds=(0.0,1.0),domain=NonNegativeReals)
     13 model.polytope2_coefficient = Var(polytope2,bounds=(0.0,1.0),domain=NonNegativeReals)
---> 14 model.polytope1 = Constraint(expr=summation(model.polytope1_coefficient) == 1, doc="Coefficients of polytope 1 must sum up to one")
     15 model.polytope2 = Constraint(expr=summation(model.polytope2_coefficient) == 1, doc="Coefficients of polytope 2 must sum up to one")
     16 model.objective = Objective(expr=max([sum([model.polytope1_coefficient[corner]*points[corner][objective] for corner in polytope1])- sum([model.polytope2_coefficient[corner]*points[corner][objective] for corner in polytope2]) for objective in range(len(points[0]))]), sense=minimize, doc='Define objective function')

/usr/lib/python2.7/site-packages/pyomo/core/base/numvalue.pyc in __eq__(self, other)
    416         (Called in response to 'self = other'.)
    417         """
--> 418         return generate_relational_expression(_eq, self, other)
    419 
    420     def __add__(self,other):

/usr/lib/python2.7/site-packages/pyomo/core/base/expr_coopr3.pyc in generate_relational_expression(etype, lhs, rhs)
   1395             etype = _eq
   1396         else:
-> 1397             raise TypeError(chainedInequalityErrorMessage())
   1398         generate_relational_expression.chainedInequality = None
   1399 

TypeError: Nonconstant relational expression used in an unexpected Boolean context.
The inequality expression:
    0.760952902864*polytope1_coefficient[0] + 0.23609925987*polytope1_coefficient[1] + 0.849584947017*polytope1_coefficient[2] - 0.00653916062779*polytope2_coefficient[3] - 0.915851403735*polytope2_coefficient[4] - 0.3996199658*polytope2_coefficient[5]  <  0.439022221545*polytope1_coefficient[0] + 0.554002150991*polytope1_coefficient[1] + 0.174480839292*polytope1_coefficient[2] - 0.204737053841*polytope2_coefficient[3] - 0.617384036536*polytope2_coefficient[4] - 0.380959902269*polytope2_coefficient[5]  <  0.752618637127*polytope1_coefficient[0] + 0.924014955365*polytope1_coefficient[1] + 0.900532369755*polytope1_coefficient[2] - 0.164859960994*polytope2_coefficient[3] - 0.826492706379*polytope2_coefficient[4] - 0.121171505813*polytope2_coefficient[5]
contains non-constant terms (variables) appearing in a Boolean context, e.g.:
    if expression <= 5:
This is generally invalid.  If you want to obtain the Boolean value of
the expression based on the current variable values, explicitly evaluate
the expression, e.g.:
    if value(expression) <= 5:
or
    if value(expression <= 5):

The fix that is suggested does not fix the problem. Instead, the code does not run at all if I add this value() function there.

Has there been any progress in one  year?

[Siirola, John D]

The problem is your use of the max() function in the objective. That construct is not representable in Pyomo expressions. You should consider adding an auxiliary variable that is constrained to be larger than each term in the max(), and then minimize that variable. 

John 

[Markus Hartikainen]

Hi,
thank you for your response. Do you mean that it does not run even once even if it does not give any errors? If I just run it once, it goes through and seems to work. What is then happening the first time?

I would have preferred to pass the min-max type problem to the solver as it is, since at least Cplex and Gurobi perform better in this kind of problems when you give them the min-max formulation instead of the linearization with the auxiliary variable. (This is based on my experiences, I do not have any real benchmarks to back this up.) Is this not possible with Pyomo?

[Markus Hartikainen]

Hi,

so I have been trying to debug this. Somehow, it seems that the problem does run the first time, because when I put tree=True, then I get the following output from glpk:

GLPSOL: GLPK LP/MIP Solver, v4.52
Parameter(s) specified in the command line:
 --write /tmp/tmpCpFlSq.glpk.raw --wglp /tmp/tmpCz3o8X.glpk.glp --cpxlp /tmp/tmpFzg7ED.pyomo.lp
Reading problem data from '/tmp/tmpFzg7ED.pyomo.lp'...
3 rows, 7 columns, 7 non-zeros
36 lines were read
Writing problem data to `/tmp/tmpCz3o8X.glpk.glp'...
35 lines were written
GLPK Simplex Optimizer, v4.52
3 rows, 7 columns, 7 non-zeros
Preprocessing...
2 rows, 4 columns, 4 non-zeros
Scaling...
 A: min|aij| =  1.000e+00  max|aij| =  1.000e+00  ratio =  1.000e+00
Problem data seem to be well scaled
Constructing initial basis...
Size of triangular part is 2
*     0: obj =  -5.671491401e-01  infeas =  0.000e+00 (0)
*     1: obj =  -7.662101680e-01  infeas =  0.000e+00 (0)
OPTIMAL LP SOLUTION FOUND
Time used:   0.0 secs
Memory used: 0.0 Mb (40420 bytes)
Writing basic solution to `/tmp/tmpCpFlSq.glpk.raw'...
12 lines were written

So, for some reason Pyomo is able to pass the problem to the solver the first time. However, it just does not work the second time around. Can somebody explain me, why this is the case.

I have been running this inside a Jupyter Notebook, so that I can check how it works. Whenever I restart the kernel, I can then solve the problem once more.

Any help is appreciated!

[Siirola, John D]

It is passing *a* problem to GLPK, but with 99.9% probability, it is not the problem you want to solve. max() is a Python function that loops over an iterable and compares each element to the running maximal value (to be completely pedantic, the initial value is compared to 0).  There (unfortunately) nothing syntactically wrong with that, so Python continues chugging along.  The only reason this is close to syntactically valid is because we allow users to express constraints as a <= b <= c, and that logic is masking the error the first time through the model.  Your comparison at the end triggers the exception hinting at the underlying problem (it is unfortunate that we cannot generate an error at the point where the original problem was raised (in the objective expression).  You should be able to see that things are off by pretty-printing your objective in your Jupyter notebook:

 

model.objective.pprint()

 

(You should see that only the last term in your max() is present).

 

Here’s the deal: Pyomo builds expressions in memory using a data structure called an expression tree, and we have no representation for “max” in the expression tree system.  Further, I am interested how you passed a bilevel (min-max) problem directly to a solver in the past?  I think you can using AMPL, but it is not clear what AMPL does with the max when it passes it into Gurobi or CPLEX (it is likely that they automatically perform the transformation for you).  I am pretty sure you can’t express something like that in an LP or MPS file…

 

john

 

 

From: pyomo...@googlegroups.com [mailto:pyomo...@googlegroups.com] On Behalf Of Markus Hartikainen
Sent: Tuesday, December 29, 2015 4:54 AM
To: Pyomo Forum
Subject: Re: [EXTERNAL] Re: Running a Concrete model multiple times in the same script

 

Hi,

--

[Markus Hartikainen]

Hi,

you are right. I tried to pretty print the objective and it is indeed just the last term in the max.

Mostly, I have been using the OPL (Optimization Programming Language) by IBM. This allows me to give problems as "min max ...". Like you said, it is also possible in ampl, but I do not know either  what it actually does when it passes that to Cplex or Gurobi. My experiences of the speed-up which I was talking about were mostly from OPL.

Having max operator is very useful especially for me, since I am often in in the situation, where I need to solve optimization problems with min max objectives. What often changes for me in these problems is that the components under the max are changing by a constant. Thus, warmstarting new problems with the optimal solutions of the previous one is really beneficial. If I use the auxiliary variables, the different problems have different feasible sets and I could not just simply use the solution as a warmstart.

However, it seems that I should add the auxiliary variable and try to do something else clever. Thank you for your fast reply!