Using a for loop for a single or specific value in variable list in constraint formulation
Answeredtime = list()
for x in range(0,int(num_t)):
time.append(x)
Resources=steel.addVars(res_list,time,vtype=GRB.INTEGER,name="Resources")
Tasks=steel.addVars(task_list,time,vtype= GRB.BINARY,name="Tasks")
My res_list is also a list
My time is also a list
Now When I am writing this,
Endresources=steel.addConstrs((gp.quicksum(Resources[res,t] for t in x)==1 for res in res_list3),name="Endresources")
where
x=[]
x.append(time[95])
But In Constraints it is showing,
Endresources[150],Endresources[151]........
Instead it should show:
Endresources[150,95],Endresources[151,95]......
How can I use "time" in Constraint instead of "x"
Endresources=steel.addConstrs((gp.quicksum(Resources[res,t] for t in x )==1 for res in res_list3),name="Endresources")
-
"The solution point you will get out of this model will tell you which constraints have to be looked at. If the value of a ArtP_constraintname or ArtN_constraintname variable is non-zero, this means that this constraint might need some relaxation or rethinking. If the value of a ArtLB_variablename or ArtUB_variablename variable is non-zero, this means that this variable might need some relaxation or rethinking of the bounds."
I haven't received any value as nonzero as you can see the feasrelax output above.Jaromil can you explain me what is happening with constraint seeing this output?
There are some weird words like:
CArtL_Resources
Can you just explain seeing those , what are the possibilities?
I have included all of those in a minimal form.
0 -
-
-
Gurobi Staff
I haven't received any value as nonzero as you can see the feasrelax output above.Jaromil can you explain me what is happening with constraint seeing this output?
You did not receive any nonzero values, because no feasible point has been found to your feasRelax model.
There are some weird words like:
Please search for CArtL_Resources in the feasRelax LP file and you will see where it appears and what it does. In partircular, CArtL_Resources is a name of an auxiliary constraint introduced by the feasRelax algorithm.
0 -
-
CArtL_Resources[1,0]: Resources[1,0] + ArtL_Resources[1,0] >= 0
This is one type of constraint :
Does that mean, there should be this type of constraint in my model which is not present?
0 -
-
-
Does that mean, there should be this type of constraint in my model which is not present?
No, this is an auxiliary constraint introduced by the feasRelax algorithm to guarantee feasibility of your model. It relaxed the constraint \(\texttt{Resources[1,0] = 0}\) to the constraint you see.
0 -
-
Hi Jaromil,
Thanks for letting me know,
I am struggling with it.
In between , I wanted to ask if you are aware of some good book on MILP programming using Python API with Gurobi.
I have already watched the MIP tutorials, apart from that are there any good books in mind?
0 -
-
Hi Jaromil,
I have one question:
What I did is : I did not changed the size of the model, Instead I disabled every constraint & ran my model my enabling them one by one.
Below are the details:
-
1. No Constraints enabled
Optimal solution found (tolerance 1.00e-04)
Best objective 0.000000000000e+00, best bound 0.000000000000e+00, gap 0.0000%-
2. End Resources Constraints enabled
Optimal solution found (tolerance 1.00e-04)
Best objective 0.000000000000e+00, best bound 0.000000000000e+00, gap 0.0000%
3. End Resources + Task execution (Process & Transfer task) constraints enabled
Optimal solution found (tolerance 1.00e-04)
Best objective 0.000000000000e+00, best bound 0.000000000000e+00,gap 0.0000% -
-
4. End Resources + Task execution (Process & Transfer task) + Energy Usage constraints enabled
Solution count 1: 18156
Optimal solution found (tolerance 1.00e-04)
Best objective 1.815600000000e+04, best bound 1.815600000000e+04, gap 0.0000%
5. End Resources + Task execution (Process & Transfer task) + Energy Usage + Task execution(group tasks) constraints enabled Solution count 1: 18870 Optimal solution found (tolerance 1.00e-04) Best objective 1.887000000000e+04, best bound 1.887000000000e+04, gap 0.0000%
Solution count 1: 18870
Optimal solution found (tolerance 1.00e-04)
Best objective 1.887000000000e+04, best bound 1.887000000000e+04, gap 0.0000%
NOTE: THERE IS ONLY ONE CONSTRAINT TO BE FILLED UP. IS THERE A PROBABLILITY OF THAT CONSTRAINT GOING WRONG RESULTING IN INFEASIBILITY?0 -
-
-
Hi Margi,
In between , I wanted to ask if you are aware of some good book on MILP programming using Python API with Gurobi.
I have already watched the MIP tutorials, apart from that are there any good books in mind?
I am not aware of any book that introduces gurobipy. To my knowledge there are only the webinars which I already commented on. You can find all recommended resources here.
IS THERE A PROBABLILITY OF THAT CONSTRAINT GOING WRONG RESULTING IN INFEASIBILITY?
I think you are on the right way. It is possible that the left out constraint is responsible for infeasibility However, it is also possible that once you add the last constraint and remove any other constraint, e.g., End Resources (or any other), that the model becomes feasible again. This would mean that the relationship between those 2 constraints may be the reason for infeasibility.
0 -
-
Hi Jaromil,
Many thanks for reply.
You are absolutely right. When I reduced my size now,(time step) & ran the program, I found that there is interconnection between Resourcebalance & Endconstraint constraint has a problem.
But I am not sure like what can be the issue ? Because I manually checked it, It seems right to me.
Below is the ouput of IIS.ilp. Can you throw some light like what can be going wrong here?
\ Model steel RTN_copy
\ LP format - for model browsing. Use MPS format to capture full model detail.
Minimize
Subject To
Endresources_time[95][141]: Resources[141,3] = 1
Resourcebalane[141,0]: Resources[141,0] = 0
Resourcebalane[141,1]: - Resources[141,0] + Resources[141,1] = 0
Resourcebalane[141,2]: - Resources[141,1] + Resources[141,2] = 0
Resourcebalane[141,3]: - Resources[141,2] + Resources[141,3] = 0
Bounds
Resources[141,0] free
Resources[141,1] free
Resources[141,2] free
Resources[141,3] free
Generals
Resources[141,0] Resources[141,1] Resources[141,2] Resources[141,3]
End0 -
-
-
Hi Margi,
Can you throw some light like what can be going wrong here?
As you probably noticed yourself, the issue is that \(\texttt{Resources[141,3]}\) is forced to 1 through the Endresources constraint and simultaneously it is forced to 0 through the chain of Resourcebalance constraints.
I think the issue lies within the Resourcebalance constraints. Your constraint make sure that the number of Resources from time n to time n+1 cannot change. But simultaneously you try to make sure that the Endresources indeed changed compared to time 0. I don't know the application of your model so I cannot really help in how to fix this.
0 -
-
Hi Jaromil,
Many thanks for reply.
I made some changes in my constraint. Also made my model size less.
It ran successfully.
Now when I am running it on my original model, it is showing me infeasible.
Will you mind letting me know for such unusual behavior?
0 -
-
-
Hi Margi,
>Will you mind letting me know for such unusual behavior?
There can be many reasons for this behavior. Maybe some constraint combination only get infeasible once you increase the time. Did you try to only slightly increasing the size of the model and checking whether it is still feasible?
Did you try now computing an IIS for the larger model?
0 -
-
Hi Jaromil,
"Did you try to only slightly increasing the size of the model and checking whether it is still feasible?"
No , I did not slightly increased, rather I made it to my original size.
"Did you try now computing an IIS for the larger model?"
I am computing it & it is still in process.
Also, How do we extract the values of variables in Gurobi after optimizing it?
0 -
-
-
You can either write the final solution to a file via
steel.write("solution.sol")or print the variable values
if steel.Status == GRB.OPTIMAL:
for v in steel.getVars():
print("%s: %f\n"%(v.VarName, v.X))0 -
-
Hi Jaromil,
First of all Huge thanks for a staff like you who responds in minimum time.
I tried my every bit in Gurobi. Lastly I got feasibility for some specific set (of data). There is an infeasibility with my original data.
Here I am posting my code. Would you please mind letting me know if the same phenomenon occurs with you too?
I have mentioned two places in bold letter where little change is required for two diff set of data.
It would really be a genuine help.
Code:
import numpy as np
import pandas as pd
import xlrd
import csv
from collections import OrderedDict
import math
import gurobipy as gp
from gurobipy import GRB
group2heats = {'1': [1, 2, 3, 4], '2': [5, 6, 7, 8], '3': [9, 10, 11, 12],
'4': [13, 14, 15, 16, 17], '5': [18, 19, 20], '6': [21, 22, 23, 24]}
equip2mw = {'EAF': 85, 'AOD': 2, 'LF': 2, 'CC1': 7, 'CC2': 7}
energy_price = [33.28, 30.00, 29.15, 28.49] #(CASE1)
#energy_price = [33.28, 30.00, 29.15, 28.49, 34.66, 50.01, #(CASE2)
# 71.52, 77.94, 81.97, 87.90, 92.76, 94.98]
# 92.31, 90.03, 90.09, 87.44, 85.37, 79.97,
# 79.92, 77.83, 76.28, 65.06, 53.07, 34.16]
equip2process_time = {
'EAF1':{'1':80.0,'2':80.0,'3':80.0,'4':80.0,'5':85.0,'6':85.0,'7':85.0,'8':85.0,'9':90.0,'10':90.0,'11':90.0,'12':90.0,'13':85.0,'14':85.0,'15':85.0,'16':85.0,'17':80.0,'18':80.0,'19':80.0,'20':80.0,'21':80.0,'22':80.0,'23':80.0,'24':80.0},
'EAF2':{'1':80.0,'2':80.0,'3':80.0,'4':80.0,'5':85.0,'6':85.0,'7':85.0,'8':85.0,'9':90.0,'10':90.0,'11':90.0,'12':90.0,'13':85.0,'14':85.0,'15':85.0,'16':85.0,'17':80.0,'18':80.0,'19':80.0,'20':80.0,'21':80.0,'22':80.0,'23':80.0,'24':80.0},
'AOD1':{'1':75.0,'2':75.0,'3':75.0,'4':75.0,'5':80.0,'6':80.0,'7':80.0,'8':80.0,'9':95.0,'10':95.0,'11':95.0,'12':95.0,'13':85.0,'14':85.0,'15':85.0,'16':85.0,'17':85.0,'18':95.0,'19':95.0,'20':95.0,'21':80.0,'22':80.0,'23':80.0,'24':80.0},
'AOD2':{'1':75.0,'2':75.0,'3':75.0,'4':75.0,'5':80.0,'6':80.0,'7':80.0,'8':80.0,'9':95.0,'10':95.0,'11':95.0,'12':95.0,'13':85.0,'14':85.0,'15':85.0,'16':85.0,'17':85.0,'18':95.0,'19':95.0,'20':95.0,'21':80.0,'22':80.0,'23':80.0,'24':80.0},
'LF1':{'1':35.0,'2':35.0,'3':35.0,'4':35.0,'5':45.0,'6':45.0,'7':20.0,'8':20.0,'9':45.0,'10':45.0,'11':45.0,'12':45.0,'13':25.0,'14':25.0,'15':25.0,'16':25.0,'17':25.0,'18':45.0,'19':45.0,'20':30.0,'21':30.0,'22':30.0,'23':30.0,'24':30.0},
'LF2':{'1':35.0,'2':35.0,'3':35.0,'4':35.0,'5':45.0,'6':45.0,'7':20.0,'8':20.0,'9':45.0,'10':45.0,'11':45.0,'12':45.0,'13':25.0,'14':25.0,'15':25.0,'16':25.0,'17':25.0,'18':45.0,'19':45.0,'20':30.0,'21':30.0,'22':30.0,'23':30.0,'24':30.0},
'CC1':{'1':50.0,'2':50.0,'3':50.0,'4':50.0,'5':60.0,'6':60.0,'7':55.0,'8':55.0,'9':60.0,'10':60.0,'11':60.0,'12':60.0,'13':70.0,'14':70.0,'15':75.0,'16':75.0,'17':75.0,'18':60.0,'19':70.0,'20':70.0,'21':50.0,'22':50.0,'23':50.0,'24':50.0},
'CC2':{'1':50.0,'2':50.0,'3':50.0,'4':50.0,'5':60.0,'6':60.0,'7':55.0,'8':55.0,'9':60.0,'10':60.0,'11':60.0,'12':60.0,'13':70.0,'14':70.0,'15':75.0,'16':75.0,'17':75.0,'18':60.0,'19':70.0,'20':70.0,'21':50.0,'22':50.0,'23':60.0,'24':60.0}
}
stage2units = {'1': {'EAF': 2}, '2': {'AOD': 2}, '3': {'LF': 2}, '4': {'CC1': 1, 'CC2': 1}}
trans_time = {'TR_S1': 10, 'TR_S2': 4, 'TR_S3': 10}
trans_time_max = {'TR_S1': 240, 'TR_S2': 240, 'TR_S3': 120}
setup_time = {'CC1': 70, 'CC2': 50}
equip2process_time['EAF'] = equip2process_time.pop('EAF1')
equip2process_time.pop('EAF2')
equip2process_time['AOD'] = equip2process_time.pop('AOD1')
equip2process_time.pop('AOD2')
equip2process_time['LF'] = equip2process_time.pop('LF1')
equip2process_time.pop('LF2')
rtn_t0 = 15
price_energy = []
for price in energy_price:
price_energy = price_energy + [int(price)] * int(60 / rtn_t0) #prices for 96 time slots
num_t = 4* 60 / rtn_t0 #total no of time slots i.e 96 #CASE 1
#num_t = 24* 60 / rtn_t0 #total no of time slots i.e 96 #CASE 2
num_stage = len((stage2units).keys()) #total no of stages i.e 4
num_groups = len((group2heats).keys()) #total no of groups i.e 6
num_heats = 0
for group, heats in group2heats.items():
num_heats += len(heats) #total no of heats i.e 24
equip2num = {} #calculating each eqipment has how many units
for unit2num in stage2units.values():
for unit, num in unit2num.items():
equip2num[unit] = num
maxPower = sum([float((equip2mw)[unit])*num for unit, num in equip2num.items()]) #(85*2 + 2*2 + 2*2 + 7*1 + 7*1)
casters = [unit for unit, num in stage2units['4'].items()]
heat_sequence = ['EAF', 'TR_EA', 'AOD', 'TR_AL', 'LF', 'TR_LC']
res_cat2idx = OrderedDict()
r_idx = 1
for stage in range(1, num_stage + 1):
for unit in (stage2units)[str(stage)].keys():
res_cat2idx[unit] = [r_idx]
r_idx = r_idx + 1
for stage in range(1, num_stage+1):
res_cat2idx['H_A_S%s' % stage] = list(range(r_idx, r_idx + num_heats)) #'H_A_S1'=EAd1 (after stage1) 'H_B_S4'=EAs4(before stage4) #'H_A_S4'=intermediate product after stage 4 i.e final products heats
r_idx = r_idx + num_heats
if int(stage) == 1:
continue
res_cat2idx['H_B_S%s' % stage] = list(range(r_idx, r_idx + num_heats))
r_idx =r_idx + num_heats
# energy resource
res_cat2idx['EN'] = [r_idx]
res_num = r_idx
[res_cat2idx, res_num]
tasks = OrderedDict()
i_idx = 1
for stage in range(1,num_stage + 1): #calculates in order how many processing and transfer tasks are there. task_num is 156 for only 3 stages.(24 heats * 6 types of tasks for 3 stages(3 processsng and 3 transfer))
for unit in (stage2units)[str(stage)].keys(): #CC1 and CC2 tasks are not counted in task_num
if int(stage) == 4:
tasks[unit] = range(i_idx, i_idx + num_groups)
i_idx = i_idx + num_groups
else:
tasks[unit] = range(i_idx, i_idx + num_heats)
i_idx = i_idx + num_heats
tasks['TR_S%s' % stage] = range(i_idx, i_idx + num_heats)
i_idx =i_idx + num_heats
task_num = i_idx - 1
[tasks, task_num]
"""calculate time slots duration of tasks"""
task_duration = dict()
task_cleanup_duration = dict()
for heat in range(1, num_heats + 1): #for heat 1, EAF,AOD,LF,TR_S1,TR_S2,TR_S3 , similar for heat 2...so on...
for task_type in ['EAF', 'AOD', 'LF']:
task_duration[tasks[task_type][heat - 1]]=int(math.ceil(float((equip2process_time)[task_type][str(heat)])/rtn_t0))
for task_type in ['TR_S1', 'TR_S2', 'TR_S3']:
task_duration[tasks[task_type][heat - 1]]=int(math.ceil(float((trans_time)[task_type])/rtn_t0))
for group in range(1, num_groups + 1):
for task_type in ['CC1', 'CC2']:
cast_time = [(equip2process_time)[task_type][str(heat)] for heat in group2heats[str(group)]]
total_time = float(sum(cast_time))
task_duration[tasks[task_type][group - 1]] = int(math.ceil(total_time/rtn_t0))
total_time = float(sum(cast_time) + setup_time[task_type])
task_cleanup_duration[tasks[task_type][group - 1]] = int(math.ceil(total_time/rtn_t0))
task_duration,task_cleanup_duration
"""calculate resource task profile, i.e. the interactions between resource and task"""
rtn_profile = dict()
for res_category, res_idxes in res_cat2idx.items(): #resources=res_cat2idx,created 174 blank dictionaries for 174 res_num
for res_idx in res_idxes:
rtn_profile[res_idx] = dict()
# equipment usage
for stage in range(1, num_stage + 1):
for unit in (stage2units)[str(stage)].keys(): #[RESOURCE][TASK]'[1][1-25]'(EAF) , '[2][49-73]'(AOD) ,'[3][97-121]'(LF),''[4][145-151]'(CC1) ,'[5][151-157](CC2)'
res = res_cat2idx[unit][0]
for task in tasks[unit]:
rtn_profile[res][task] = [-1] + [0] * (task_duration[task] - 1) + [1]
if stage == 4:
rtn_profile[res][task] = [-1] + [0] * (task_cleanup_duration[task] - 1) + [1]
# heat consumption and generation for the first three stages
for heat_idx in range(0, num_heats):
# process task generate intermediate heat H_A_S (after stage)
for task_cat, res_cat in [('EAF', 'H_A_S1'), ('AOD', 'H_A_S2'), ('LF', 'H_A_S3')]:
task = tasks[task_cat][heat_idx]
resource = res_cat2idx[res_cat][heat_idx]
rtn_profile[resource][task] = [0] * task_duration[task] + [1]
# process task consume intermediate heat H_B_S (before stage)
for task_cat, res_cat in [('AOD', 'H_B_S2'), ('LF', 'H_B_S3')]:
task = tasks[task_cat][heat_idx]
resource = res_cat2idx[res_cat][heat_idx]
rtn_profile[resource][task] = [-1] + [0] * task_duration[task]
# transfer task transports heat
for task_cat, res_cat in [('TR_S1', ['H_A_S1', 'H_B_S2']),('TR_S2', ['H_A_S2', 'H_B_S3']),('TR_S3', ['H_A_S3', 'H_B_S4'])]:
task = tasks[task_cat][heat_idx]
resource1 = res_cat2idx[res_cat[0]][heat_idx]
rtn_profile[resource1][task] = [-1] + [0] * task_duration[task]
resource2 =res_cat2idx[res_cat[1]][heat_idx]
rtn_profile[resource2][task] = [0] * task_duration[task] + [1]
# group-heat consumption and generation
heat_consume_time_in_group = dict()
heat_generate_time_in_group = dict()
for unit in (stage2units)['4']:
heat_consume_time_in_group[unit] = dict()
heat_generate_time_in_group[unit] = dict()
for group, heats in group2heats.items():
task = tasks[unit][int(group)-1]
duration = 0
for heat in heats:
consume_time = int(math.floor(duration/rtn_t0))
resource = res_cat2idx['H_B_S4'][heat - 1]
rtn_profile[resource][task] = [0] * (task_duration[task] + 1)
rtn_profile[resource][task][consume_time] = -1
duration += (equip2process_time)[unit][str(heat)]
generate_time = int(math.ceil(duration/rtn_t0))
resource = res_cat2idx['H_A_S4'][heat - 1]
rtn_profile[resource][task] = [0] * (task_duration[task] + 1)
rtn_profile[resource][task][generate_time] = 1
heat_consume_time_in_group[unit][heat] = consume_time
heat_generate_time_in_group[unit][heat] = generate_time
# energy usage
for stage, units in stage2units.items():
stage = int(stage)
for unit in units.keys():
norm_mw = float(equip2mw[unit])
for task in tasks[unit]:
total_energy = 0
if stage < 4:
heat = task - tasks[unit][0] + 1
total_energy = norm_mw * equip2process_time[unit][str(heat)] / 60
elif stage == 4:
group = task - tasks[unit][0] + 1
for heat in group2heats[str(group)]:
total_energy += norm_mw * equip2process_time[unit][str(heat)] / 60
profile = [norm_mw * rtn_t0 / 60] * (task_duration[task] - 1)
profile = profile + [total_energy - sum(profile)] + [0]
rtn_profile[res_cat2idx['EN'][0]][task] = profile
rtn_profile, heat_consume_time_in_group, heat_generate_time_in_group
time = list()
for x in range(0,int(num_t)):
time.append(x)
res_list = [res for cat, cat_res in res_cat2idx.items() for res in cat_res]
task_list=[task for cat, cat_task in tasks.items() for task in cat_task]
steel=gp.Model('steel RTN')
Resources=steel.addVars(res_list,time,vtype=GRB.INTEGER,name="Resources")
Tasks=steel.addVars(task_list,time,vtype= GRB.BINARY,name="Tasks")
heat2tasks = [task for task_cat, task_list2 in tasks.items()
if 'CC' not in task_cat for task in task_list2]
Taskexecution=steel.addConstrs((gp.quicksum(Tasks[task,t] for t in time)==1 for task in heat2tasks),name="Taskexecution")
y={}
for res_cat, resources in res_cat2idx.items():
if 'H_A_S' not in res_cat or 'H_A_S4' in res_cat:
continue
y[res_cat]= list(resources)
for key in y:
y_list = y[key]
Transfertime=steel.addConstrs((Resources[res,t] == 0 for res in y_list for t in time),name ="Transfertime")
res_list3 = [res for cat, cat_res in res_cat2idx.items() if cat == 'H_A_S4' for res in cat_res]
Endresources=steel.addConstrs((Resources[res,time[(int(num_t)-1)]] ==1 for res in res_list3 ), name="Endresources_time[95]")
group2tasks = dict()
for group in range(1, num_groups + 1):
group2tasks[group] = [tasks[caster][group - 1] for caster in stage2units['4'].keys()]
for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for t in time ) ==1
for task in group2tasks[key]),name="Taskexecution1")
x={}
for res_cat, resources in res_cat2idx.items():
if 'H_B_' not in res_cat:
continue
max_wait_time = trans_time_max['TR_S%d' % (int(res_cat[-1]) - 1)]
min_tran_time = trans_time['TR_S%d' % (int(res_cat[-1]) - 1)]
for res in resources: # each heat a constraint
x[res] = (math.ceil((max_wait_time-min_tran_time)/rtn_t0))
Transfertime1= steel.addConstrs((gp.quicksum(Resources[res,t] for t in time )<=x[res] for res in x.keys()),name="Transfertime1")
res_list1 = [res for cat, cat_res in res_cat2idx.items() if cat != 'EN' for res in cat_res]
Resourcebalance = {}
for res in res_list1:
for t in time:
if t == 0:
if 1 <= res and res <= 3:
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] == 2, name="Resourcebalane[%s,%d]"%(res,t))
elif 4<= res and res<=5:
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] == 1, name="Resourcebalane[%s,%d]"%(res,t))
#else:
# Resourcebalance[res,t] = steel.addConstr(Resources[res,t] == 0, name="Resourcebalane[%s,%d]"%(res,t))
else:
sumexpression = gp.LinExpr(0)
for task in task_list:
if task not in rtn_profile[res]:
continue
duration = len(rtn_profile[res][task]) - 1
for theta in range(0, min(duration, t) + 1):
sumexpression.add(rtn_profile[res][task][theta]* Tasks[task,t-theta] )
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] == Resources[res,time[time.index(t)-1]] + sumexpression, name="Resourcebalane[%s,%d]"%(res,t))
res = res_cat2idx['EN'][0]
Energyusage = {}
for ti in range(int(num_t)):
sumexpression = gp.LinExpr(0)
for task_cat, task_list1 in tasks.items():
if 'TR' in task_cat: # transportation has no energy consumption
continue
for task in task_list1:
duration = task_duration[task]
for theta in range(0, min(duration, ti) + 1):
sumexpression.add(rtn_profile[res][task][theta]* Tasks[task,ti-theta] )
Energyusage[res,ti] = steel.addConstr(Resources[res,ti] ==sumexpression,name="Energyusage[%s,%d]"%(res,ti))
reslist2=[]
reslist2.append(res_cat2idx['EN'][0])
obj=gp.quicksum(price_energy[t]*Resources[res,t] for res in reslist2 for t in time)
steel.setObjective(obj,GRB.MINIMIZE)
steel.optimize()I tried to get every possible information but I think now I am on the correct path but missing some little bit.#
Huge thanks for the help.
0 -
-
-
It looks like your initial values for \(\texttt{Resources[4,0], Resources[5,0]}\) are too strict. Relaxing the constraints to inequalities makes the model feasible
[...]
elif 4<= res and res<=5:
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] >= 1, name="Resourcebalance[%s,%d]"%(res,t))Please note that I don't know whether this is feasible for your application. It is just one way to make your model feasible.
With the No Relaxation heuristic, a feasible point is found quite quickly. You can use the heuristic by setting the corresponding parameter NoRelHeurTime. 200 seconds was enough on my machine.
steel.setParam("NoRelHeurTime", 200)
steel.optimize()Note that the model is quite complex and may take a lot of time to reach the default MIPGap of 0.01%. You might want to increase the MIPGap if your application allows for this. The model might also solve just fine with the default MIPGap. Alternatively, you can just set a TimeLimit.
0 -
-
Hi Jaromil,
Many thanks for the reply.
1. Did you ran the model on case 2?
2. Also, there is one constraint associated with those two resources(Res 4 & Res 5). I am listing that below, Let me know if I have formulated in right manner.
Constraint explanation & formulation:
The formulation:
group2tasks = dict()
for group in range(1, num_groups + 1):
group2tasks[group] = [tasks[caster][group - 1] for caster in stage2units['4'].keys()]
for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for t in time ) ==1
for task in group2tasks[key]),name="Taskexecution1")3. "With the No Relaxation heuristic, a feasible point is found quite quickly. You can use the heuristic by setting the corresponding parameter NoRelHeurTime. 200 seconds was enough on my machine"
- What is no relaxation heuristic? 200 sec was the total running time of model?
I kept my model to run for almost 8 hours- It still did not got completed.
4.Why would model take long time to reach MIPGap?
- Huge thanks Jaromil.
Kind regards,
Margi
0 -
-
-
Hi Margi,
1. Did you ran the model on case 2?
I took the code you posted adjusted the Resourcebalance constraints as mentioned in my previous comment and ran it for
energy_price = [33.28, 30.00, 29.15, 28.49]
[...]
num_t = 4 * 60 / rtn_t0energy_price = [33.28, 30.00, 29.15, 28.49, 34.66, 50.01,
71.52, 77.94, 81.97, 87.90, 92.76, 94.98]
[...]
num_t = 4 * 60 / rtn_t0and
energy_price = [33.28, 30.00, 29.15, 28.49, 34.66, 50.01,
71.52, 77.94, 81.97, 87.90, 92.76, 94.98
92.31, 90.03, 90.09, 87.44, 85.37, 79.97,
79.92, 77.83, 76.28, 65.06, 53.07, 34.16]
[...]
num_t = 24* 60 / rtn_t0The first 2 cases are definitely feasible. The first 2 solved in under 600 seconds with the NoRelHeurTime=200 parameter set. The last, bigger case takes way longer, which makes sense, because it is 6x larger. It did not solve yet but I think that it is feasible as well. You might need to increase the NoRelHeurTime parameter by a lot for this one, maybe 3600 is enough.
Also, there is one constraint associated with those two resources(Res 4 & Res 5). I am listing that below, Let me know if I have formulated in right manner.
It looks like the constraint is missing a sum. Currently your constraint says
\[\begin{align*}
\sum_{t \in time} Tasks_{task,t} = 1 \quad\forall task \in group2tasks
\end{align*}\]
So the \(\sum_{u\in CCs}\) is missing. The Tasks variable does not even have a \(u\) index.What is no relaxation heuristic? 200 sec was the total running time of model?
It is a heuristic which aims to find a feasible solution without ever computing a root node relaxation. You have to explicitly set it.
steel.setParam("NoRelHeurTime",200)
steel.optimize()This will run the heuristic for 200 seconds.
I kept my model to run for almost 8 hours- It still did not got completed.
Did you set the no relaxation heuristic parameter? Could you show a log snippet? The first ~30 lines and last ~20 lines are enough.
Why would model take long time to reach MIPGap?
When a model is complex, finding feasible points and proving a strong dual bound can get very hard. As described in the description of the MIPGap parameter, the MIPGap defines the relative difference between the best feasible point and a proven best bound. It can be hard to reach a very low MIPGap for large complex models. You might want to have a look at our MIP tutorials to learn more about basics of MIPs.
Best regards,
Jaromił0 -
-
-
Hi Margi,
The largest model, i.e., the one with
energy_price = [33.28, 30.00, 29.15, 28.49, 34.66, 50.01,
71.52, 77.94, 81.97, 87.90, 92.76, 94.98
92.31, 90.03, 90.09, 87.44, 85.37, 79.97,
79.92, 77.83, 76.28, 65.06, 53.07, 34.16]
[...]
num_t = 24* 60 / rtn_t0is infeasible in its current state. However, it can be made feasible by making the Resources variables continuous instead of discrete.
Resources=steel.addVars(res_list,time,vtype=GRB.CONTINUOUS,name="Resources")
on top of the additional change of
if 1 <= res and res <= 3:
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] == 2, name="Resourcebalance[%s,%d]"%(res,t))
elif 4<= res and res<=5:
Resourcebalance[res,t] = steel.addConstr(Resources[res,t] >= 1, name="Resourcebalance[%s,%d]"%(res,t))You don't have to use any special parameters in this case.
I don't know whether making the Resource variables continuous makes sense for your application so it's up to you to understand what is going wrong.
0 -
-
Hi Jaromil,
First of all huge thanks for the reply.
Let me go step by step.
1. "It looks like the constraint is missing a sum"
Yes, I think the same because ,
The code looks like:
group2tasks = dict()
for group in range(1, num_groups + 1):
group2tasks[group] = [tasks[caster][group - 1] for caster in stage2units['4'].keys()]
for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for t in time ) ==1
for task in group2tasks[key]),name="Taskexecution1")
where: output:
group2tasks:
{1: [145, 151],
2: [146, 152],
3: [147, 153],
4: [148, 154],
5: [149, 155],
6: [150, 156]}
Taskexecution1[145]: Tasks[145,0] + Tasks[145,1] + Tasks[145,2]
+ Tasks[145,3] + Tasks[145,4] + Tasks[145,5] + ...
Taskexecution1[151]: Tasks[151,0] + Tasks[151,1] + Tasks[151,2]
+ Tasks[151,3] + Tasks[151,4] +...
Instead it should look like:
Tasks[145,1]+Tasks[151,1]=1
Tasks[146,1]+Tasks[152,1]=1
2)Did you set the no relaxation heuristic parameter? Could you show a log snippet? The first ~30 lines and last ~20 lines are enough.
No Jaromil, I did not put Heuristic.
3)I don't know whether making the Resource variables continuous makes sense for your application so it's up to you to understand what is going wrong.
I think I can try that but initially,
let me put the below changes :
- Adding Hueristics
-Adding MIPGap
-Reformulation of constraint(Task execution)
Please let me know what are your views on this?
Also Is it possible to reformulate that constraint in the manner I mentioned?
- Kind Regards
Margi0 -
-
-
I think I can try that but initially,
let me put the below changes :
- Adding Hueristics
-Adding MIPGap
-Reformulation of constraint(Task execution)Experimenting with the heuristic or the MIPGap parameter are only optional and you should not waste too much time on this. You should focus on formulating your model correctly first. After you take care of the Taskexecution1 constraint, you should go through every constraint 1 by 1 and check whether they are correct before proceeding.
Instead it should look like:
Tasks[145,1]+Tasks[151,1]=1
Tasks[146,1]+Tasks[152,1]=1Given the formula you have shown, I don't think it's correct, because you still have only one sum \(\sum_{task \in group2tasks} Tasks_{task,1} = 1 \), i.e., you are doing something different with the time as is done in the formula.
The above constraint can be formulated viafor key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for task in group2tasks[key] ) ==1 for t in time),
name="Taskexecution1[%d]"%key)You can use the above snippet to adjust the constraint in a way that you think it is correct. If you want to only use \(t=1\) as in your example, you can just say \(\texttt{for t in [1]}\) instead of \(\texttt{for t in time}\). I think from here you should be able to proceed on your own.
0 -
-
Many thanks Jaromil for the reply.
"After you take care of the Taskexecution1 constraint, you should go through every constraint 1 by 1 and check whether they are correct before proceeding"
I think you are completely right. I will now recheck my every constraint very carefully one by one.
Before that I have one simple que related to Gurobi:
Does Gurobi support conditional constraints?
By that what I mean is:
group2tasks:
{1: [145, 151],
2: [146, 152],
3: [147, 153],
4: [148, 154],
5: [149, 155],
6: [150, 156]}for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for t in time ==1 for task in group2tasks[value1 or value 2] ),
name="Taskexecution1[%d]"%key)Can I use such conditions to get the below output?
For value1 of group2tasks[1]:
Tasks[145,0]+Tasks[145,1]+......+Tasks[145,95] = 1 OR
For value2 of group2tasks[1]:
Tasks[151,0]+Tasks[151,1]+......+Tasks[151,95] = 1Does Gurobi support such conditional constraints from given set of data?
0 -
-
This constraint seems to be very typical constraint.
I searched necessary documentation on Gurobi website but still how can I evaluate any one specific constraint from given set of data?
I tried to do via the for loop,
i.e for group,tasks in group2tasks.items():
Jaromil, Is it Ideally possible doing so?
0 -
-
-
Community Moderator
Hi Margi,
it is possible to model conditional statements. For more information, have a look at Jaromił's entry in knowledge base.
Best regards
Jonasz0 -
-
Hi,
Thanks for the response.
Jaromil has helped me in learning Gurobi as this is completely new to me. Its Overwhelming.
However I think I have a question which I am unable to find on any of the resources of Gurobi.
"Any of the one constraint from the given set of two constraint should be executed."
By this , What I mean:
group2tasks:
{1: [145, 151],
2: [146, 152],
3: [147, 153],
4: [148, 154],
5: [149, 155],
6: [150, 156]}for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for task in group2tasks[key] ) ==1 for t in time),
name="Taskexecution1[%d]"%key)Output is:
Taskexecution1[1][0] = Tasks[145,0]+Tasks[151,0] =1
Taskexecution1[1][0] = Tasks[145,1]+Tasks[151,1] =1.....
What it should be:
Taskexecution1[1] = Tasks[145,0]+Tasks[145,1] + Tasks[145,2]+Tasks[145,3] +...=1
OR
Taskexecution1[1] = Tasks[151,0]+Tasks[151,1] + Tasks[151,2]+Tasks[151,3] +...=1Jaromil, Jaromił Najman
Would you please kindly let mw know the direction to proceed on such issue?
I would really appreciate it.
I am almost there. I have manually checked my each and every constraint and it is perfectly right.
0 -
-
for key in group2tasks:
Taskexecution1=steel.addConstrs(
(gp.quicksum(Tasks[task,t] for t in time ) ==1
for task in group2tasks[key][0] or
for task in group2tasks[key][1]),name="Taskexecution1")I have tried to formulate in this manner.
0 -
-
b = steel.addVar(vtype=GRB.BINARY, name="b")
for key in group2tasks:
Taskexecution1= steel.addConstr((b == 1) >> ((gp.quicksum(Tasks[task,t] for t in time ) ==1
for task in group2tasks[key])), name="indicator_constr1")
for key in group2task:
Taskexecution2= steel.addConstr((b == 0) >> ((gp.quicksum(Tasks[task,t] for t in time ) ==1 for task in group2tasks[key][1])), name="indicator_constr1")Jaromił Najman , Have I used Indictor Constraints in a correct manner? I have tried to do this too, still I lack in execution of this constraint.
0 -
-
-
Hi Margi,
Your indicator constraint usage is almost correct. It should be
for key in group2tasks:
b = steel.addVar(vtype=GRB.BINARY, name="b[%d]"%key)
Taskexecution1= steel.addConstr((b == 1) >> (gp.quicksum(Tasks[group2tasks[key][0],t] for t in time ) ==1),
name="indicator_constr1[%d]"%key)
Taskexecution2= steel.addConstr((b == 0) >> (gp.quicksum(Tasks[group2tasks[key][1],t] for t in time ) ==1),
name="indicator_constr2[%d]"%key)Note multiple things. I am using only one \(\texttt{for}\)-loop to have access to the same \(\texttt{key}\) value when introducing the two indicator constraints. I am also introducing a new binary variable for every indicator constraint pair. This is necessary if you want to have an OR decision for every \(\texttt{key}\) value. If you really want to have only one binary variable control all OR decisions at once, then you have to put the addition of the auxiliary binary variable \(b\) before the \(\texttt{for}\)-loop. I am also not using an additional
for task in group2tasks[key]
in the addConstr method. This is just not necessary and makes the code more complicated and harder to read. You know which element of the \(\texttt{group2tasks[key]}\) list you want to access so you can just do it directly. I used the \(\texttt{key}\) value when naming the binary variable and constraints to avoid duplicate names. It is very important to avoid duplicate names when working with an LP file.
0 -
-
Hi Jaromil,
First of all, thanks. I think now with your help I am slowly learning how to use Gurobi.
However, in matter of Indicator Constraint I have one doubt:
Your post "How do I model conditional statements in Gurobi? – Gurobi Help Center" states that:
"Gurobi's indicator constraints require a binary variable as the indicator variable, we model if x>y by enforcing x>y→b=1 and x≤y→b=0. The binary variable b thus indicates if x>y is true (b=1) or false (b=0)"
I am not sure if I require to model my constraint as an Indicator constraint. I am stating down below my formulation, let me know if it requires an indicator assignment as ["x>y→b=1"] there are no x & y variables in my model.
Below is my formulation:
group2tasks
{1: [145, 151], 2: [146, 152], 3: [147, 153], 4: [148, 154], 5: [149, 155], 6: [150, 156]}
Taskexecution1 = {}
for key in group2tasks:
y_list = group2tasks[key]
Taskexecution1[key]=steel.addConstr((gp.quicksum(Tasks[task,t] for t in time for task in y_list) ==1),name="Taskexecution1[%d]"%key)This is the output I am getting:
Taskexecution1[1]: Tasks[145,0] + Tasks[145,1] + Tasks[145,2]
+ Tasks[145,3] + Tasks[145,4] + Tasks[145,5] + Tasks[145,6]
+ Tasks[145,7] + Tasks[145,8] + Tasks[145,9] + Tasks[145,10]
+ Tasks[145,11] + Tasks[145,12] + Tasks[145,13] + Tasks[145,14]
+ Tasks[145,15] + Tasks[151,0] + Tasks[151,1] + Tasks[151,2]
+ Tasks[151,3] + Tasks[151,4] + Tasks[151,5] + Tasks[151,6]
+ Tasks[151,7] + Tasks[151,8] + Tasks[151,9] + Tasks[151,10]
+ Tasks[151,11] + Tasks[151,12] + Tasks[151,13] + Tasks[151,14]
+ Tasks[151,15] = 1However, what I am expecting is
[Tasks[145,0]+.....+Tasks[145,t] ] + [Tasks[151,0]+.....+Tasks[151,t] ] == 1
i.e Either my first constraint (or part of constraint) ([Tasks[145,0]+.....+Tasks[145,t] ] ) should be 1
OR
Second part of constraint [Tasks[151,0]+.....+Tasks[151,t] ] should be 1
i.e From 1st key of group2tasks, either 145 should be executed or 151 should be executed.
I am not sure whether I should construct it with Indicator constraint method as I did before or like adding two binary variable constraint in Gurobi?
I want to know If my understanding goes in a right manner or not.
Kind regards,
Margi
0 -
-
-
The Knowledge Base article How do I model conditional statements in Gurobi? – Gurobi Help Center is meant to help in modeling \(\texttt{if}\)-conditions with \(x > y\). You don't have such a condition, thus you don't have to follow the Knowledge Base article.
The usage of indicator constraints as formulated in my previous comment should be enough in your case. What the code
for key in group2tasks:
b = steel.addVar(vtype=GRB.BINARY, name="b[%d]"%key)
Taskexecution1= steel.addConstr((b == 1) >> (gp.quicksum(Tasks[group2tasks[key][0],t] for t in time ) ==1),
name="indicator_constr1[%d]"%key)
Taskexecution2= steel.addConstr((b == 0) >> (gp.quicksum(Tasks[group2tasks[key][1],t] for t in time ) ==1),
name="indicator_constr2[%d]"%key)models is the following
\[\begin{align*}
\text{if } b = 1 \text{ then } \sum_t Tasks[145,t] =1\\
\text{if } b = 0 \text{ then } \sum_t Tasks[151,t] =1
\end{align*}\]Since, \(b\) is a binary variable, it will be either \(0\) or \(1\) in the final solution. Thus, only one of the above constraints will be active in the final solution. In other words, either 145 or 151 will be executed. Note that the above does not mean that if \(b=1\) then \(\sum_t Tasks[151,t] \neq 1\). If \(b=1\) then \(\sum_t Tasks[151,t]\) can take any value, which could also be \(1\). To me it looks like you would like to force \(\sum_t Tasks[151,t] = 0 \) when \(b=1\). You can enforce this, with 2 different constraints.
for key in group2tasks:
b = steel.addVar(vtype=GRB.BINARY, name="b[%d]"%key)
Taskexecution1= steel.addConstr((gp.quicksum(Tasks[group2tasks[key][0],t] for t in time ) == b),
name="constr1[%d]"%key)
Taskexecution2= steel.addConstr((gp.quicksum(Tasks[group2tasks[key][1],t] for t in time ) == 1 - b),
name="constr2[%d]"%key)The above models
\[\begin{align*}
\sum_t Tasks[145,t] &= b\\
\sum_t Tasks[151,t] &= 1-b
\end{align*}\]So you see that when \(b=1\) then \(\sum_t Tasks[151,t] = 0\) and if \(b=0\) then \(\sum_t Tasks[151,t] = 1\).
A different question: This post and this post looks extremely similar to your latest comment. Did you post it? If yes, could you please delete it to avoid duplicates?
0 -
-
Hi Jaromil,
Many thanks for making me clear regarding the concept of Indicator Constraints. I have few doubts with regards to it if you can answer:
First of all, briefing out the constraint formulation & the expected output:
Below if what should be expected from it:
For group I : Tasks[145,0]+...+Tasks[145,t] ==1 OR Tasks[151,0]+...+Tasks[151,t] ==1
So, the question is,
1. Why would we require TWO constraints to execute this single constraint?
2. Also, if we are adding a new binary variable, in the model who will decide if b=0 or b=1?
3. For such type of constraint, is it compulsory to use binary variable?
Many thanks for such a clear explanation.
And yes, the posts are similar. I will delete them.
0 -
Post is closed for comments.
Comments
210 comments