Groups keyboard shortcuts have been updated
Dismiss
See shortcuts

how can I put variables into if statement using python, pyomo and a linear solver(Gurobi)?

240 views
Skip to first unread message

Armaghan Bhr

unread,
Jul 18, 2020, 12:34:17 PM7/18/20
to Pyomo Forum
0

I'm implementing an optimization problem in pyomo. However, I need to define a constraint based on temperatures variables by using if statement. The constraint is like: if self.b.temperatures['supply',t]=self.b.temperatures['return',t] then: self.b.temperatures['supply',t]=self.b.temperatures['supply',t-1] else:
self.b.heat_flow[t]=(self.b.temperatures['supply',t]-self.b.temperatures['return',t])*self.cp * b.mass_flow[t]

I know that I cannot use a variable in if statement. So, now I'm looking for some tricks to still implement this constraint. A part of the code is as below:

self.block.heat_flow = Var(self.TIME, within=NonNegativeReals) lines = self.params['lines'].v()

        def _mass_flow(b, t):
           
return self.params['mass_flow'].v(t)

       
self.block.mass_flow = Param(self.TIME, rule=_mass_flow)

       
       
def _decl_init_heat_flow(b):
           
return b.heat_flow[0] == (
                   
self.params['temperature_supply'].v() -
                   
self.params['temperature_return'].v()) * \
                   
self.cp * b.mass_flow[0]

       
self.block.decl_init_heat_flow = Constraint(
            rule
=_decl_init_heat_flow)

       
self.block.temperatures = Var(lines, self.TIME)




Daniel Gebbran

unread,
Jul 19, 2020, 8:08:23 PM7/19/20
to Pyomo Forum
Hello Armaghan,

I have recently struggled through a similar problem as you.
In short, you are correct: you cannot declare an "if" constraint based on the value of the variable. Pyomo constructs the model based on values prior to executing the solver, so you need expressions that have been defined instead of any expressions based on variables. You can, however, try some things, depending on the structure of the problem.

First, you could try to implement a piecewise function (there are straight up implementations of piecewise functions within Pyomo, there's documentation on their website and some examples!) if it fits your model (maybe it doesn't, just by looking quickly at what you've written down). 
Second, you can try to analyze the values of the variables and re-declare the constraint (which we can further discuss if it fits your model) if you could solve your process in iterations, instead of solving in a centralized manner.

Best,
Daniel Gebbran

Armaghan Bhr

unread,
Jul 20, 2020, 6:32:05 AM7/20/20
to Pyomo Forum
Hi Daniel,

Thanks a lot for your explanation!

Regarding your solutions, I'm familiar with piecewise functions in pyomo but I only used it to linearize some product terms. Moreover, I cannot use the second method( re-declaring the constraints) as it's an optimization problem and I solve it at each time step(no iteration).

Well, I've tried to implement big-M expression to cover if statement. Also, I considered some more constraints to make the product of binary and continuous variable linear as I'm using a linear solver. However, I don't get a reasonable result. Would you please check if it's correct or at least close to what you've done before?


self.block.del_T_on = Var(self.TIME, within=Binary)

# linearization of the product term b.del_T_on[t]*b.temperatures['supply', t-1] by defining inequality constraints:
self.block.temp_del_on = Var(self.TIME)
def _temp_on_1_u(b, t):
   
return b.temp_del_on[t] <= self.params['temperature_max'].v() * b.del_T_on[t]

def _temp_on_1_l(b, t):
   
return b.temp_del_on[t] >= self.params['temperature_min'].v() * b.del_T_on[t]

def _temp_on_2_u(b, t):
   
return b.temp_del_on[t] <= self.params['temperature_max'].v()

def _temp_on_2_l(b, t):
   
return b.temp_del_on[t] >= min(self.params['temperature_min'].v(), 0)

def _temp_on_3_l(b, t):
   
if t==0:
       
return Constraint.Skip
   
else:
       
return b.temp_del_on[t] >= b.temperatures['supply', t-1] - (1 - b.del_T_on[t]) * self.params['temperature_max'].v()

def _temp_on_3_u(b, t):
   
if t==0:
       
return Constraint.Skip
   
else:
       
return b.temp_del_on[t] <= b.temperatures['supply', t-1] - (1 - b.del_T_on[t]) * self.params['temperature_min'].v()

def _temp_on_4(b, t):
   
if t==0:
       
return Constraint.Skip
   
else:
       
return b.temp_del_on[t] <= b.temperatures['supply', t-1] + (1 - b.del_T_on[t]) * self.params['temperature_max'].v()

self.block.temp_on_1_u = Constraint(self.TIME, rule=_temp_on_1_u)
self.block.temp_on_1_l = Constraint(self.TIME, rule=_temp_on_1_l)
self.block.temp_on_2_u = Constraint(self.TIME, rule=_temp_on_2_u)
self.block.temp_on_2_l = Constraint(self.TIME, rule=_temp_on_2_l)
self.block.temp_on_3_u = Constraint(self.TIME, rule=_temp_on_3_u)
self.block.temp_on_3_l = Constraint(self.TIME, rule=_temp_on_3_l)
self.block.temp_on_4 = Constraint(self.TIME, rule=_temp_on_4)


# linearization of the product term b.del_T_on[t]*b.temperatures['return', t] by defining inequality constraints:
self.block.temp_del_on_r = Var(self.TIME)
def _temp_on_r_1_u(b, t):
   
return b.temp_del_on_r[t] <= self.params['temperature_return_max'].v() * b.del_T_on[t]

def _temp_on_r_1_l(b, t):
   
return b.temp_del_on_r[t] >= self.params['temperature_return_min'].v() * b.del_T_on[t]

def _temp_on_r_2_u(b, t):
   
return b.temp_del_on_r[t] <= self.params['temperature_return_max'].v()

def _temp_on_r_2_l(b, t):
   
return b.temp_del_on_r[t] >= min(self.params['temperature_return_min'].v(), 0)

def _temp_on_r_3_l(b, t):
   
return b.temp_del_on_r[t] >= b.temperatures['return', t] - (1 - b.del_T_on[t]) * \
           
self.params['temperature_return_max'].v()

def _temp_on_r_3_u(b, t):
   
return b.temp_del_on_r[t] <= b.temperatures['return', t] - (1 - b.del_T_on[t]) * \
           
self.params['temperature_return_min'].v()

def _temp_on_r_4(b, t):
   
return b.temp_del_on[t] <= b.temperatures['return', t] + (1 - b.del_T_on[t]) * self.params['temperature_return_max'].v()

self.block.temp_on_r_1_u = Constraint(self.TIME, rule=_temp_on_r_1_u)
self.block.temp_on_r_1_l = Constraint(self.TIME, rule=_temp_on_r_1_l)
self.block.temp_on_r_2_u = Constraint(self.TIME, rule=_temp_on_r_2_u)
self.block.temp_on_r_2_l = Constraint(self.TIME, rule=_temp_on_r_2_l)
self.block.temp_on_r_3_u = Constraint(self.TIME, rule=_temp_on_r_3_u)
self.block.temp_on_r_3_l = Constraint(self.TIME, rule=_temp_on_r_3_l)
self.block.temp_on_r_4 = Constraint(self.TIME, rule=_temp_on_r_4)

# linearization of the product term b.del_T_on[t]*b.heat_flow[t] by defining inequality constraints:
self.block.heat_del_on = Var(self.TIME)
def _heat_on_1(b, t):
   
return b.heat_del_on[t] <= b.Qmax * b.del_T_on[t]

def _heat_on_2(b, t):
   
return b.heat_del_on[t] <= b.heat_flow[t]

def _heat_on_3(b, t):
   
return b.heat_del_on[t] >= b.heat_flow[t] - (1 - b.del_T_on[t]) * b.Qmax

def _heat_4(b, t):
   
return 0 <= b.heat_flow[t]

self.block.heat_on_1 = Constraint(self.TIME, rule=_heat_on_1)
self.block.heat_on_2 = Constraint(self.TIME, rule=_heat_on_2)
self.block.heat_on_3 = Constraint(self.TIME, rule=_heat_on_3)
self.block.heat_4 = Constraint(self.TIME, rule=_heat_4)


###################################################################

def _decl_temperatures(b, t):
   
if t == 0:
       
return Constraint.Skip
   
elif b.mass_flow[t] == 0:
     
return Constraint.Skip
   
else:
       
return b.temperatures['supply', t] - (b.temperatures['supply', t - 1] - b.temp_del_on[
            t
]) == b.heat_del_on[t] / b.mass_flow[t] / self.cp + b.temp_del_on_r[t]

self.block.decl_temperatures = Constraint(self.TIME, rule=_decl_temperatures)


Daniel Bacilla

unread,
Jul 21, 2020, 7:35:49 PM7/21/20
to pyomo...@googlegroups.com
Hi Armaghan,

The thing is, specifying a constraint within an if statement is doable as long as the statement is not based on a variable's value.
Pyomo will build the constraint based on the statement at the time it builds the model... and that means it may read the variable's value (which may be 0, 1, or any other value upon initialization, or just trash) and build the constraint based on that value, at the moment it initializes the model. 
Declaring a piecewise linear constraint (as well as using a BigM model) may circumvent that... but it then requires the usage of binary/integer variables. If you then circumvent the usage of binary variables with yet additional constraints, you're reducing even further the success ratio for your solution. What are the chances you can use a mixed-integer solver? 



--
You received this message because you are subscribed to a topic in the Google Groups "Pyomo Forum" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/pyomo-forum/ZujyGkozTQk/unsubscribe.
To unsubscribe from this group and all its topics, send an email to pyomo-forum...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/pyomo-forum/14a21062-b086-4607-96d8-8d8cbb757b97o%40googlegroups.com.

Armaghan Bhr

unread,
Jul 22, 2020, 6:44:51 PM7/22/20
to Pyomo Forum
Hi Daniel,

I get your point. I use Gurobi solver which covers both linear problems and mixed-integer linear problems.

Well, here I have binary variable self.block.del_T_on and I tried to implement "  if self.b.temperatures['supply',t]=self.b.temperatures['return',t] then: self.b.temperatures['supply',t]=self.b.temperatures['supply',t-1] else:
self.b.heat_flow[t]=(self.b.temperatures['supply',t]-self.b.temperatures['return',t])*self.cp * b.mass_flow[t]" in a form of one constraint (self.block.decl_temperatures). However, as the slover is linear and I had the product term of binary and continuous variables in that equation, I needed to define more constraints on variables self.block.temp_del_on, self.block.temp_del_on_r and self.block.heat_del_on but it seems the if statement isn't implemented correctly, yet.

To unsubscribe from this group and all its topics, send an email to pyomo...@googlegroups.com.

Daniel Bacilla

unread,
Jul 22, 2020, 8:37:37 PM7/22/20
to pyomo...@googlegroups.com
Hello Armaghan,

Just so I understand this right, both  self.b.temperatures['supply',t]  and   self.b.temperatures['return',t]  are variables, correct?
The problem is (IMO) not the following product of binary and continuous variables, but the "if" statement using values of variables, as I previously discussed.   

To unsubscribe from this group and all its topics, send an email to pyomo-forum...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/pyomo-forum/932e8e08-0833-4390-81b8-b695201e1aeao%40googlegroups.com.

Armaghan Bhr

unread,
Jul 23, 2020, 12:18:55 AM7/23/20
to Pyomo Forum
Hi Daniel,

Yes. self.b.temperatures['supply',t]  and self.b.temperatures['return',t]  are continuous variables. 

Well, I've somehow combined the 2 equations: "self.b.temperatures['supply',t]=self.b.temperatures['supply',t-1]" and
"self.b.heat_flow[t]=(self.b.temperatures['supply',t]-self.b.temperatures['return',t])*self.cp * b.mass_flow[t]". So, first, it was like: " b.temperatures['supply', t] - (1-b.del_T_on[t]) * b.temperatures['supply', t - 1] == (b.del_T_on[t]) * (b.heat_flow[t] / b.mass_flow[t] / self.cp + b.temperatures['return', t])" and I defined those 3 variables(self.block.temp_del_on, self.block.temp_del_on_r and self.block.heat_del_on) and their constraints to be substituted for the product terms of "b.del_T_on[t] * b.temperatures['supply', t - 1]", "b.del_T_on[t] * (b.heat_flow[t]" and "b.del_T_on[t] * b.temperatures['return', t]". Hence, if the binary var is 1, then "self.b.heat_flow[t]=(self.b.temperatures['supply',t]-self.b.temperatures['return',t])*self.cp * b.mass_flow[t]" is satisfied and if the binary var=0, I'll get "self.b.temperatures['supply',t]=self.b.temperatures['supply',t-1]". However, my problem is that still sth is missing in defining the if statement in my code and I don't know what other eq/factor I need to take into account.

Daniel Bacilla

unread,
Jul 23, 2020, 7:25:19 AM7/23/20
to pyomo...@googlegroups.com
Hi Armaghan,

Thanks for confirming their use as variables.
In short, the "something is missing" you are looking for is, again: you cannot declare a conditional (if) constraint based on the values of something (the variables) that has not yet been assigned.

Think of it like that: you want to pay for bread, but you don't tell the baker how many breads you are taking home. He will now know how many breads to give you, not how much he should charge you for. You need to say a value X for those breads. In this case, the number X is the values for both of your variables self.b.temperatures['supply',t]  and self.b.temperatures['return',t].

And likewise, you want the constraint to be defined based on an if statement, but it does not know how many breads it will have.
Because, after all, self.b.temperatures['supply',t]  and self.b.temperatures['return',t] are results of your optimization problem. So, there is no point on your if constraint.
Even if we can express these conditions mathematically, they cannot be implemented computationally in this kind of solver. And interpreters like Pyomo cannot build this constraint, simply put.
I say this from my own experience :) I have a problem in which I would like to have, similarly, a conditional constraint based on a variable value yet to be calculated... which cannot be implemented - I'm still trying to somehow solve it. If I eventually can, I'll let you know!

Finally, let me leave you with this:  https://groups.google.com/forum/#!topic/pyomo-forum/SJpBPHZlk7Q which is exactly what you are asking for. So, there are ways around it, and (unfortunately) I could not yet implement it, but others have done so - have a look at this other post. 
A last piece of advice is to use: value(self.b.temperatures['supply',t]), when you are reading the actual value for any variable or constraint outside the model per se. Not sure if that helps, but it surely won't hurt if you didn't know about it yet.

Hope I could help!
Best regards,
Daniel Gebbran

To unsubscribe from this group and all its topics, send an email to pyomo-forum...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/pyomo-forum/cf52e8cb-e873-48e4-9707-f359f7f2d522o%40googlegroups.com.

Armaghan Bhr

unread,
Aug 1, 2020, 3:50:00 PM8/1/20
to Pyomo Forum
Thanks a lot, Daniel for clarifying!

It really gave me a clear idea of this subject.

Best regards,
Armaghan

Soheil.mt

unread,
Aug 2, 2020, 1:52:50 AM8/2/20
to pyomo...@googlegroups.com
Hi Armaghan,

i had same problem when i was modeled a piecewise-linear function in pyomo as below:
image.png
where:
image.png
so, i defined a new binary variable image.png and redefined constrain as bellow:

image.png
where M is a big number.
and coded in pyomo:

model.W1 = Var(model.V1, within = PositiveReals)              ### Arrival Time at 1st-Level Nudes Variabel
model.W2 = Var(model.V2, within = PositiveReals)              ### Arrival Time at 2ed-Level Nudes Variabel
model.e = Param(model.Vc)                                     ### Time Window Lower Bond
model.l = Param(model.Vc)                                     ### Time Window Upper Bond
model.alpha = Param(model.Vc)                                 ### Prefer Time
model.P = Var(model.VcT, within = Binary)
        
def Fx (model, j):
    W = model.W1[j] if j in model.Vc3 else model.W2[j]
    if model.e[j]!= model.alpha[j]:
        return (W - model.e[j])/(model.alpha[j] - model.e[j])
    else:
        return 0
        
def Gx (model, j):
    W = model.W1[j] if j in model.Vc3 else model.W2[j]
    if model.l[j]!= model.alpha[j]:
        return (model.l[j] - W)/(model.l[j] - model.alpha[j])
    else:
        return 0
        
def __Expr1(model, j):
    W = model.W1[j] if j in model.Vc3 else model.W2[j]
    return W - model.alpha[j] <= M * (1-model.P[j])
model.Expr1 = Constraint(model.VcT , rule = __Expr1)
        
def __Expr2(model, j):
    W = model.W1[j] if j in model.Vc3 else model.W2[j]
    return model.alpha[j] - W <= M * model.P[j]
model.Expr2 = Constraint(model.VcT , rule = __Expr2)
        
def __Expr3(model, j):
    return Fx(model, j) * model.P[j] + Gx(model, j) * (1-model.P[j])
model.Mu = Expression(model.VcT , rule = __Expr3)


I hope this will help.





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.
To view this discussion on the web visit https://groups.google.com/d/msgid/pyomo-forum/3c5641a7-1588-4b29-a372-ad1745a8c3a5o%40googlegroups.com.

Ahmad Heidari

unread,
Feb 26, 2025, 7:44:29 PMFeb 26
to Pyomo Forum
Hello Armaghan,

I hope you are doing well.

Suppose X, Y, and Z are your variables; each of which could be functions of other sets. For example, X[i,t], Y[i,t], and Z[i,t].
And, suppose we have the following statement: 

if      X[i,t] == Xmax[i]    and    Y[i,t] > Ymax[i]     then     Z[i,t] == Y[i,t] - Ymax[i] 

To solve this problem, Pyomo uses Generalized Disjunctive Programs (GDPs) as follows.

# Small epsilon
epsilon = 1e-3

# Rule IF
def if_disjunct_rule(disj, I, T):
    disj.con1 = Constraint(expr = MP0.X[I,T] - MP0.Xmax[I]  <= epsilon)
    disj.con2 = Constraint(expr = -MP0.X[I,T] + MP0.Xmax[I] <= epsilon)
    disj.con3 = Constraint(expr = MP0.Y[I,T] - MP0.Ymax[I]  >= epsilon)
    disj.con4 = Constraint(expr = MP0.Z[I,T] == MP0.Y[I,T] - MP0.Ymax[I])

# Rule Else
def else_disjunct_rule(disj, I, T):
    disj.inactive_on = Constraint(expr = MP0.Z[I,T] == 0)

# Create disjuncts for the Rules above
MP0.I4_if      = Disjunct(MP0.I, MP0.T, rule=if_disjunct_rule)
MP0.I4_else = Disjunct(MP0.I, MP0.T, rule=else_disjunct_rule)

# For each (I,T), one of the two disjuncts must hold.
def disjunction_rule(MP0, I, T):
    return [MP0.I4_if[I,T], MP0.I4_else[I,T]]
MP0.I4_disjunction = Disjunction(MP0.I, MP0.T, rule=I4_disjunction_rule)

Finally, before your solve statement, please provide the code below:

# The bigM value depends on your problem!
TransformationFactory('gdp.bigm').apply_to(MP0, bigM=1000000)
 solver = SolverFactory('gurobi')

For further reading, you may refer to Pyomo Book on their website: https://www.pyomo.org/documentation
or, you can refer to the tutorial book that we've published: https://scholarsmine.mst.edu/gradstudent_works/6/
Reply all
Reply to author
Forward
0 new messages