Using a for loop which has range of values in the constraint formulation
AnsweredI am trying to formulate below constraint.
I have modified it as:
I have calculated the LHS part as per some data in my code which is as below:
res_cat2idx= [OrderedDict([('EAF', [1]),
('AOD', [2]),
('LF', [3]),
('CC1', [4]),
('CC2', [5]),
('H_A_S1', range(6, 30)),
('H_A_S2', range(30, 54)),
('H_B_S2', range(54, 78)),
('H_A_S3', range(78, 102)),
('H_B_S3', range(102, 126)),
('H_A_S4', range(126, 150)),
('H_B_S4', range(150, 174)),
('EN', [174])])
trans_time={'TR_S1': 10, 'TR_S2': 4, 'TR_S3': 10}
trans_time_max = {'TR_S1': 240, 'TR_S2': 240, 'TR_S3': 120}
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.append(math.ceil((max_wait_time-min_tran_time)/rtn_t0))
Now I am trying to formulate it this way, where [x] is my LHS part. However I am getting below error:
Transfertime= steel.addConstrs((gp.quicksum(Resources[res,t] for t in time )==x for res in resources),name="Transfertime")
Kindly help me with this,
Regards,
Margi
-
Hi Jaromił Najman,
Thank you so much for your help so far. Your comments were really useful. I identified the error, it was due to some another constraint.
My doubt is :
Suppose for my objective function:obj=gp.quicksum(price_energy[t-1]*EnergyResource[res,t] for res in e_res for t in time1) - gp.quicksum(sp_price[t-1]*SpinningResource[res,t] for res in s_res for t in time1)
The 1st term:
gp.quicksum(price_energy[t-1]*EnergyResource[res,t] for res in e_res for t in time1)
represents the EN cost.
and the 2nd term:gp.quicksum(sp_price[t-1]*SpinningResource[res,t] for res in s_res for t in time1)
represents the SP revenue
My question is : Is it possible in Gurobi to know the values of these terms in the objective function after the optimisation has been performed?
Just like below:
Thank you.
Margi.0 -
-
-
Gurobi Staff
Hi Margi,
Good to hear that you are making progress.
My question is : Is it possible in Gurobi to know the values of these terms in the objective function after the optimisation has been performed?
Yes, there are two ways to get those values. Both ways are performed after the optimization process was successful.
First one:
You can use solution information and recompute the values:
if model.status == GRB.OPTIMAL:
EN_value = 0
SP_value = 0
for t in time1:
for res in e_res:
EN_value += price_energy[t-1]*EnergyResource[res,t].X
for res in s_res:
SP_value +=sp_price[t-1]*SpinningResource[res,t].X
print(f"EN_value: {EN_value}")
print(f"SP_value: {SP_value}")Second one:
You can add auxiliary variable to catch these values and let Gurobi do the job
aux_EN_val = model.addVar(lb=-GRB.INFINITY, name="aux_EN_val")
model.addConstr(aux_EN_val == gp.quicksum(price_energy[t-1]*EnergyResource[res,t] for res in e_res for t in time1))
aux_SP_val = model.addVar(lb=-GRB.INFINITY, name="aux_SP_val")
model.addConstr(aux_SP_val == gp.quicksum(sp_price[t-1]*SpinningResource[res,t] for res in s_res for t in time1))
...
model.optimize()
if model.status == GRB.OPTIMAL:
print(f"EN_value: {aux_EN_val.X}")
print(f"SP_value: {aux_SP_val.X}")Here there is the downside of introducing redundant variables and equalities to the model, but as long as this is not done too often, it is fine and should not have any impact on solver's performance.
For convenience reasons I would prefer the second option. However, both options are correct and it is up to you to choose one.
Best regards,
Jaromił0 -
-
Hi Jaromił Najman,
Thank you so much. I'll try to implement the second approach in my code.
0 -
-
Hi Jaromił Najman,
Suppose these are price values for 96 time slots:price_energy = [33, 33, 33, 33, 30, 30, 30, 30, 29, 29, 29, 29, 28, 28, 28, 28, 34, 34, 34, 34, 50, 50, 50, 50, 71, 71, 71, 71, 77, 77, 77, 77, 81, 81, 81, 81, 87, 87, 87, 87, 92, 92, 92, 92, 94, 94, 94, 94, 92, 92, 92, 92, 90, 90, 90, 90, 90, 90, 90, 90, 87, 87, 87, 87, 85, 85, 85, 85, 79, 79, 79, 79, 79, 79, 79, 79, 77, 77, 77, 77, 76, 76, 76, 76, 65, 65, 65, 65, 53, 53, 53, 53, 34, 34, 34, 34]
& Suppose these are the generation values for 96 time slots:
E_Wind=[7.623491061,4.517699791,3.555646106,2.48438853,0.391585596,1.421687974,1.860133449,0.686599093,
2.214638619,2.942800387,3.988794377,6.969082667,9.23995314,10.64126725,10.46401467,10.11317679,
10.25090409,9.420465543,8.118983344,8.649518668,8.125910457,8.615698059,10.46931187,14.98334437,
19.50308155,17.39479448,14.5269699,12.61508684,14.02984771,15.76651556,17.33734019,17.38501503,
19.98105231,21.98991494,27.85554933,35.53486477,36.97855651,37.74094637,40,39.9535476,39.66668364,
37.53965263,29.76783986,23.27020832,17.78026792,14.89940406,11.82376611,10.21626853,9.290887791,
7.289359751,5.651708858,6.812203942,8.832883411,10.47705394,12.23491061,11.11801559,11.76182957,
14.3696837,12.55396526,11.41751133,11.80461468,13.22385779,12.81026843,12.90847043,10.34095655,
9.97952427,8.952681709,10.87271431,8.957571436,8.596954108,8.731014109,9.608720012,11.57846483,
12.4961035,13.05271736,11.58457699,10.10910202,10.91916671,11.02225844,10.63556257,8.29908827,
7.467427291,6.326491112,5.687566852,4.413793103,6.339530383,5.93286813,2.206081597,1.344674782,2.710130902,
2.766362757,2.951764886,3.475373096,1.560230225,0.514236235,0]
Now what I do is generate 50 random numbers for price as well as generation for each time slot.import numpy as np
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]
# Generate 50 random numbers for each mean in energy_price with variance 5
random_numbers_dict = {}
for mean in energy_price:
random_numbers = np.random.normal(loc=mean, scale=np.sqrt(5), size=50)
random_numbers_dict[mean] = random_numbers.tolist()
# Print the dictionary
for mean, random_numbers in random_numbers_dict.items():
print(f"Mean: {mean}, Random Numbers: {random_numbers}")So I intend to perform optimisation at each sample for each mean.
So I reframed this constraint as:
Before it was:for res in res_list2:
for t in time1:
# add constr aux = EnergyResouce - E
steel.addConstr(aux[res,t] == EnergyResource[res,t] - E_Wind[t-1],name="aux_constr_%s_%d"%(res,t))
# add constr maxaux = max(0, aux)
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))Now it is :
for mean, random_numbers in random_numbers_dict.items():
for number in random_numbers:
for res in res_list2:
for t in time1:
# add constr aux = EnergyResouce - E
steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d"%(res,t))
# add constr maxaux = max(0, aux)
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))The final optimisation loop is :
for mean, random_numbers in random_numbers_dict.items():
for number in random_numbers:
obj=gp.quicksum(number*maxaux[res,t] for res in e_res for t in time1)
steel.setObjective(obj,GRB.MINIMIZE)
steel.optimize()Does it seems correct in terms of Gurobi Format?
Also how do I get objective function value at each mean, Meaning :-
At each price --> At each random value ---> There is some objective value.
So at each price --> There will be 50 objective function values, How do I retrive and store it at each iteration of price ? Because I require it to perform further calculation?
I hope it is clear. If not, please let me know.
Your guidance will be very important.
Best regards,
Margi0 -
-
-
Hi Margi,
I think it should be
for mean in random_numbers_dict:
reandom_numbers = random_numbers_dict[mean] for number in random_numbers:You can always test Python only stuff yourself by just using \(\texttt{print}\) statements and running the code, e.g.,
for mean, random_numbers in random_numbers_dict.items(): for number in random_numbers:
print(number)and see whether you get all numbers you want.
Also how do I get objective function value at each mean, Meaning :-
At each price --> At each random value ---> There is some objective value.
So at each price --> There will be 50 objective function values, How do I retrive and store it at each iteration of price ? Because I require it to perform further calculation?After the optimization process is finished, you can retrieve the objective value by accessing the ObjVal attribute.
steel.setObjective(obj,GRB.MINIMIZE) steel.optimize()
objval = -GRB.INFINITY
if (steel.status == GRB.OPTIMAL):
objval = steel.ObjVal
else:
print("CANNOT RETRIEVE OBJVAL")
# do something with the objval valueBest regards,
Jaromił0 -
-
Hi Jaromił Najman,
Thank you for your input.
The changes in the model I have is just addition of for loop.
The code is below:price_energy = [33, 33, 33, 33, 30, 30, 30, 30, 29, 29, 29, 29, 28, 28, 28, 28, 34, 34, 34, 34, 50, 50, 50, 50, 71, 71, 71, 71, 77, 77, 77, 77, 81, 81, 81, 81, 87, 87, 87, 87, 92, 92, 92, 92, 94, 94, 94, 94, 92, 92, 92, 92, 90, 90, 90, 90, 90, 90, 90, 90, 87, 87, 87, 87, 85, 85, 85, 85, 79, 79, 79, 79, 79, 79, 79, 79, 77, 77, 77, 77, 76, 76, 76, 76, 65, 65, 65, 65, 53, 53, 53, 53, 34, 34, 34, 34]
E_Wind=[7.623491061,4.517699791,3.555646106,2.48438853,0.391585596,1.421687974,1.860133449,0.686599093,
2.214638619,2.942800387,3.988794377,6.969082667,9.23995314,10.64126725,10.46401467,10.11317679,
10.25090409,9.420465543,8.118983344,8.649518668,8.125910457,8.615698059,10.46931187,14.98334437,
19.50308155,17.39479448,14.5269699,12.61508684,14.02984771,15.76651556,17.33734019,17.38501503,
19.98105231,21.98991494,27.85554933,35.53486477,36.97855651,37.74094637,40,39.9535476,39.66668364,
37.53965263,29.76783986,23.27020832,17.78026792,14.89940406,11.82376611,10.21626853,9.290887791,
7.289359751,5.651708858,6.812203942,8.832883411,10.47705394,12.23491061,11.11801559,11.76182957,
14.3696837,12.55396526,11.41751133,11.80461468,13.22385779,12.81026843,12.90847043,10.34095655,
9.97952427,8.952681709,10.87271431,8.957571436,8.596954108,8.731014109,9.608720012,11.57846483,
12.4961035,13.05271736,11.58457699,10.10910202,10.91916671,11.02225844,10.63556257,8.29908827,
7.467427291,6.326491112,5.687566852,4.413793103,6.339530383,5.93286813,2.206081597,1.344674782,2.710130902,
2.766362757,2.951764886,3.475373096,1.560230225,0.514236235,0]
# Generate 50 random numbers for each mean in price_energy with variance 5
electricity_random_numbers_dict = {}
for mean in price_energy:
random_numbers = np.random.normal(loc=mean, scale=np.sqrt(5), size=50)
electricity_random_numbers_dict[mean] = random_numbers.tolist()
# Generate 50 random numbers for each mean in E_Wind with variance 5
generation_random_numbers_dict = {}
for mean in E_Wind:
shape = 2
# Scale parameter calculation
scale = mean / (stats.gamma(1 + 1/shape).mean()**(1/shape))
random_numbers = np.random.weibull(shape, size=50) * scale
generation_random_numbers_dict[mean] = random_numbers.tolist()Constraint:
for mean in generation_random_numbers_dict:
random_numbers = generation_random_numbers_dict[mean]
for number in random_numbers:
for res in res_list2:
for t in time1:
# add constr aux = EnergyResouce - E
steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d"%(res,t))
# add constr maxaux = max(0, aux)
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))Objective function:
Objective_values =[]
for mean in electricity_random_numbers_dict:
random_numbers =electricity_random_numbers_dict[mean]
for number in random_numbers:
obj=gp.quicksum(number*maxaux[res,t] for res in e_res for t in time1)
steel.setObjective(obj,GRB.MINIMIZE)
steel.optimize()
objval = -GRB.INFINITY
if (steel.status == GRB.OPTIMAL):
objval = steel.ObjVal
Objective_values.append(objval)
else:
print("CANNOT RETRIEVE OBJVAL")Result for each optimisation run:
Gurobi Optimizer version 11.0.0 build v11.0.0rc2 (win64 - Windows 10.0 (19045.2)) CPU model: 12th Gen Intel(R) Core(TM) i7-1255U, instruction set [SSE2|AVX|AVX2]count: 10 physical cores, 12 logical processors, using up to 12 threads Optimize a model with 477504 rows, 23739 columns and 1032788 nonzeros Model fingerprint: 0xf736a10b Model has 460800 general constraints Variable types: 288 continuous, 23451 integer (11520 binary) Coefficient statistics: Matrix range [2e-01, 2e+01] Objective range [3e+01, 3e+01] Bounds range [1e+00, 2e+00] RHS range [4e-02, 8e+01] Presolve removed 0 rows and 126 columns Presolve time: 11.96s Explored 0 nodes (0 simplex iterations) in 12.01 seconds (2.99 work units) Thread count was 1 (of 12 available processors) Solution count 0 Model is infeasible Best objective -, best bound -, gap - CANNOT RETRIEVE OBJVAL Gurobi Optimizer version 11.0.0 build v11.0.0rc2 (win64 - Windows 10.0 (19045.2))
This repeats for every optimisation run.
Did I represent anything incorrect? Can you throw some light why model suddenly became infeasible?Best Regards,
Margi0 -
-
-
Your model is infeasible. You need to find out why. Please refer to How do I determine why my model is infeasible? in this case.
0 -
-
Hi Jaromił Najman,
Thank you for directing me. I understand to check the constraint by writing steel.write("steel.lp").I followed the same.
I am recognizing the problem with the below constraint:
Before it was:
for res in res_list2:
for t in time1:
# add constr aux = EnergyResouce - E
steel.addConstr(aux[res,t] == EnergyResource[res,t] - E_Wind[t-1],name="aux_constr_%s_%d"%(res,t))
# add constr maxaux = max(0, aux)
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))After it is :
for res in res_list2:
for t in time1: #loop for time slots
for mean in generation_random_numbers_dict: # loop for each generation value
random_numbers = generation_random_numbers_dict[mean] # getting 50 random numbers
for number in random_numbers:
print(number)
steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d"%(res,t))
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))I checked,I mentioned the loop correctly.
Initially for each t, there was only one generation value.
Now for each t, there are 50 generation values, and I loop them in the innermost loop.
I double checked the .lp files.
I think aux_constr_124_1 is repeating 50 times. How do I take care of it within the loop ?Is there a way in Gurobi to take care of this naming conventions? because I think I mentioned the constraint correctly , but the issue is because of the naming conventions?
Best regards,
Margi0 -
-
-
Is it possible that you put the loop
for mean in generation_random_numbers_dict: # loop for each generation value random_numbers = generation_random_numbers_dict[mean] # getting 50 random numbers for number in random_numbers:
over every constraint set that you build and additionally over the objective and optimize statement?
If yes, then this means that you add all constraints to one and the same model and then optimize over it with different objectives multiple times. You should rather implement something like the following pseudo code
def buildModel(number):
# your code which constructs the whole model and optimizes it
...
objval = -GRB.INFINITY
if (model.status == GRB.OPTIMAL): objval = model.ObjVal
return objval else: print("CANNOT RETRIEVE OBJVAL")
exit()
# start of the big loop
for mean in generation_random_numbers_dict: # loop for each generation value
random_numbers = generation_random_numbers_dict[mean] # getting 50 random numbers
for number in random_numbers:
print(f"Running model for mean {mean} and number {number}")
objval = buildModel(number)
Objective_values.append(objval)0 -
-
Hi Jaromił Najman,
All other constraints in my model doesn't require 'number' in their formulation, even though do they require to put the loop around them?
For ex:y={}
for res_cat, resources in res_cat2idx.items():
if 'H_A_S' not in res_cat or 'H_A_S3' in res_cat:
continue
y[res_cat]= list(resources)
Transfertime = {}
for mean in generation_random_numbers_dict:
random_numbers = generation_random_numbers_dict[mean]
for number in random_numbers:
for key in y:
y_list = y[key]
Transfertime=steel.addConstrs((Resources[res,t] == 0 for res in y_list for t in time1),name ="Transfertime")The only constraint that requires 'number' in the formulation is below:
for res in res_list2: for t in time1: #loop for time slots for mean in generation_random_numbers_dict: # loop for each generation value random_numbers = generation_random_numbers_dict[mean] # getting 50 random numbers for number in random_numbers: print(number) steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d"%(res,t)) steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))
Regards,
Margi
0 -
-
-
Hi Margi,
All other constraints in my model doesn't require 'number' in their formulation, even though do they require to put the loop around them?
Here you have the freedom of choice.
You can put the loops around everything and always build the model from scratch. This is the easy method and I would recommend it as long as your model building process does not take too long.
You can add all unaffected constraints to the model only exactly once. You can then put the loops only around the affected constraints and objective and adjust them in every iteration. In this case, you again have two options. You will either have to remove all affected constraints after the optimization process, because you will add new ones in the next loop iteration. Or you add the constraints without the number value once at the beginning, and then in the loops, you only adjust the RHS attribute of each affected constraint. For that you should save the old number value so you can first subtract/add it to the RHS and then add/subtract the new number value.
Best regards,
Jaromił0 -
-
Hi Jaromił Najman,
Thank you. I'll go with your recommended option.
The only doubt is regarding the affected constraint. I am sure, I mentioned the constraint correctly. But I think instead of 1 constraint in each time slot, now there are 50 with the same names. It is true, as per logic.But due to same naming conventions, there is a warning message. Will it affect the feasibility?
Constraint:for res in res_list2:
for t in time1:
time_mean_dict[t] = mean
random_numbers = generation_random_numbers_dict[mean]
for number in random_numbers:
steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d"%(res,t))
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t)).LP file:
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0
aux_constr_124_1: - EnergyResource[124,1] + aux[124,1] = 0It repeats for maxaux constraint also. So each of those is repeated 50 times.
Warning:
Is there a way to handle it ?
0 -
-
-
The only doubt is regarding the affected constraint. I am sure, I mentioned the constraint correctly. But I think instead of 1 constraint in each time slot, now there are 50 with the same names. It is true, as per logic.
Yes, there will be 50 constraints with the same name. Are you sure that you want to add all 50? In the output you posted it looks like \(\texttt{number=0}\) for all cases.
But due to same naming conventions, there is a warning message. Will it affect the feasibility?
No, this will not affect feasibility. It will just make retrieving constraint information harder.
Is there a way to handle it ?
Sure, as long as \(\texttt{number}\) is unique, you can use this in the name.
for number in random_numbers: steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d_%f"%(res,t, number)) steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d_%f"%(res,t, number))
If \(\texttt{number}\) is not unique, you can add a counter
counter = 0
for number in random_numbers:
steel.addConstr(aux[res,t] == EnergyResource[res,t] - number ,name="aux_constr_%s_%d_%d_%f"%(res,t, counter,mean))
steel.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d_%d_%f"%(res,t, counter,mean))
counter += 10 -
-
Hi Jaromił Najman,
I am facing a very weird error. I want to understand where I am going wrong.
So I checked the constraint formulation of SpinningResource. It is what it exactly should be.
The cost/objective value I am getting through Gurobi is 76. I retrieved the solution and tried to calculate the cost by hand and according to constraint formulation also it should be 73. I am not sure why Gurobi is not reaching the optimal. The optimality gap is 0.01%.
For ex: According to constraint formulation of Gurobi: The value of SpinningResource at 85th time step should be 34. But when I see the solution file, It is 0.
Can you guide me what I could check to understand this error?
0 -
-
-
The cost/objective value I am getting through Gurobi is 76. I retrieved the solution and tried to calculate the cost by hand and according to constraint formulation also it should be 73. I am not sure why Gurobi is not reaching the optimal. The optimality gap is 0.01%.
Could you please post the log output of the run?
How exactly do you retrieve the objective value and how do you calculate it by hand?
You could write the model to an LP file and share it. Note that uploading files in the Community Forum is not possible but we discuss an alternative in Posting to the Community Forum.
Best regards,
Jaromił0 -
-
Hi Jaromił Najman,
Please find the log file :
Gurobi Optimizer version 11.0.0 build v11.0.0rc2 (win64 - Windows 10.0 (19045.2)) CPU model: 12th Gen Intel(R) Core(TM) i7-1255U, instruction set [SSE2|AVX|AVX2] Thread count: 10 physical cores, 12 logical processors, using up to 12 threads Optimize a model with 17042 rows, 23837 columns and 125422 nonzeros Model fingerprint: 0xb5dd7353 Model has 96 general constraints Variable types: 386 continuous, 23451 integer (11520 binary) Coefficient statistics: Matrix range [2e-01, 9e+01] Objective range [2e+00, 9e+01] Bounds range [1e+00, 2e+00] RHS range [4e-01, 4e+01] Presolve removed 12550 rows and 13547 columns Presolve time: 0.58s Presolved: 4492 rows, 10290 columns, 78880 nonzeros Variable types: 94 continuous, 10196 integer (6112 binary) Root simplex log... Iteration Objective Primal Inf. Dual Inf. Time 28170 7.5378400e+04 3.876600e+02 0.000000e+00 5s 28521 7.5388408e+04 0.000000e+00 0.000000e+00 5s Root relaxation: objective 7.538841e+04, 28521 iterations, 4.45 seconds (6.27 work units) Nodes | Current Node | Objective Bounds | Work Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time 0 0 75388.4080 0 925 - 75388.4080 - - 5s H 0 0 118184.73649 75388.4080 36.2% - 7s 0 0 75931.9942 0 864 118184.736 75931.9942 35.8% - 9s H 0 0 95200.736489 75932.3238 20.2% - 9s H 0 0 95199.736489 75932.3238 20.2% - 9s 0 0 75958.8549 0 781 95199.7365 75958.8549 20.2% - 10s 0 0 75964.4529 0 874 95199.7365 75964.4529 20.2% - 11s 0 0 75964.7616 0 894 95199.7365 75964.7616 20.2% - 11s 0 0 75964.7616 0 905 95199.7365 75964.7616 20.2% - 11s 0 0 75964.7648 0 886 95199.7365 75964.7648 20.2% - 11s H 0 0 93981.653155 75964.7686 19.2% - 12s 0 0 75964.7686 0 855 93981.6532 75964.7686 19.2% - 12s 0 0 75965.5734 0 895 93981.6532 75965.5734 19.2% - 12s 0 0 75965.6137 0 900 93981.6532 75965.6137 19.2% - 12s 0 0 75965.6148 0 914 93981.6532 75965.6148 19.2% - 12s H 0 0 93965.819822 75965.6148 19.2% - 13s 0 0 75965.6148 0 893 93965.8198 75965.6148 19.2% - 13s 0 0 75966.0380 0 946 93965.8198 75966.0380 19.2% - 13s 0 0 75966.0380 0 829 93965.8198 75966.0380 19.2% - 13s 0 2 75966.0380 0 824 93965.8198 75966.0380 19.2% - 14s 1 4 75966.0380 1 813 93965.8198 75966.0380 19.2% 229 15s H 35 46 92058.967966 75967.3872 17.5% 620 18s H 36 46 76936.555536 75967.3872 1.26% 603 18s H 37 46 76372.705929 75967.3872 0.53% 587 18s H 81 86 76370.705929 75967.3872 0.53% 445 19s 86 96 75970.1743 10 992 76370.7059 75967.3872 0.53% 439 20s 159 169 75972.8183 14 839 76370.7059 75967.3872 0.53% 319 26s H 160 169 76370.705856 75967.3872 0.53% 317 26s H 163 169 76369.372596 75967.3872 0.53% 315 26s 286 303 75991.5663 19 894 76369.3726 75967.3872 0.53% 230 30s 367 317 75994.3067 21 977 76369.3726 75967.3872 0.53% 292 37s 378 330 75994.5230 21 937 76369.3726 75967.3872 0.53% 362 43s 414 359 76124.4064 23 998 76369.3726 75967.3872 0.53% 463 47s 492 415 cutoff 25 76369.3726 75967.3872 0.53% 507 50s 543 442 cutoff 24 76369.3726 75967.3872 0.53% 580 56s 570 442 cutoff 34 76369.3726 75969.3184 0.52% 608 62s H 573 442 76369.372512 75969.3184 0.52% 609 62s 587 458 76137.2492 6 894 76369.3725 75969.3184 0.52% 638 65s 656 478 76351.9009 8 889 76369.3725 75969.3184 0.52% 684 72s 712 487 cutoff 11 76369.3725 75969.7346 0.52% 706 75s 799 547 76123.0049 9 785 76369.3725 75971.1676 0.52% 766 82s 863 564 76135.3730 7 985 76369.3725 75971.1676 0.52% 801 86s 948 621 76143.0581 10 708 76369.3725 75971.1676 0.52% 843 92s 1017 663 cutoff 12 76369.3725 75971.1676 0.52% 858 97s 1096 728 76206.7599 12 645 76369.3725 75971.1676 0.52% 866 101s 1200 809 76238.1776 13 640 76369.3725 75971.1676 0.52% 846 105s 1440 920 76368.5010 25 484 76369.3725 75974.0394 0.52% 798 112s 1586 998 76287.1487 8 701 76369.3725 75974.0394 0.52% 766 116s 1737 1030 cutoff 9 76369.3725 75974.0394 0.52% 740 120s 1810 1064 cutoff 9 76369.3725 75974.0394 0.52% 740 125s 1883 1114 cutoff 9 76369.3725 75974.0394 0.52% 757 132s 2001 1126 cutoff 27 76369.3725 75978.0542 0.51% 763 144s 2038 1141 76230.9392 20 910 76369.3725 75978.0542 0.51% 796 153s 2088 1209 76326.4799 21 821 76369.3725 75978.0542 0.51% 829 163s 2198 1230 cutoff 20 76369.3725 75979.4320 0.51% 840 171s 2261 1260 cutoff 15 76369.3725 75979.4320 0.51% 865 181s 2339 1266 76207.9812 15 843 76369.3725 75979.4320 0.51% 878 191s H 2352 1266 76369.372493 75979.4320 0.51% 888 191s 2366 1367 76302.0489 17 754 76369.3725 75979.4320 0.51% 894 201s 2535 1368 76208.2649 12 829 76369.3725 75979.4320 0.51% 877 205s 2538 1370 76347.2978 31 755 76369.3725 75979.4320 0.51% 876 212s 2544 1374 76367.3247 44 762 76369.3725 75993.8835 0.49% 873 216s 2558 1383 76362.3598 17 813 76369.3725 76000.5854 0.48% 869 220s 2565 1388 76040.0682 6 773 76369.3725 76032.4202 0.44% 866 225s 2577 1396 76363.7521 41 814 76369.3725 76044.6414 0.43% 862 230s 2584 1401 76198.5720 27 943 76369.3725 76045.4803 0.42% 860 235s 2592 1406 76354.3949 22 917 76369.3725 76048.4102 0.42% 857 241s 2601 1413 76357.9323 50 929 76369.3725 76048.6455 0.42% 885 247s 2602 1414 76365.7715 55 738 76369.3725 76048.6522 0.42% 885 251s 2608 1418 76352.3371 20 946 76369.3725 76048.6522 0.42% 883 255s 2619 1425 76353.0790 43 836 76369.3725 76054.9104 0.41% 879 261s 2625 1429 76354.5087 21 801 76369.3725 76055.6989 0.41% 877 265s 2633 1435 76283.9236 14 918 76369.3725 76055.6989 0.41% 874 270s 2658 1460 76057.0014 28 969 76369.3725 76056.8965 0.41% 896 275s 2688 1478 76063.4682 29 880 76369.3725 76056.9583 0.41% 899 280s 2709 1481 76063.7322 30 877 76369.3725 76056.9583 0.41% 903 285s 2751 1467 76121.1501 32 843 76369.3725 76056.9583 0.41% 906 291s 2803 1436 76170.2600 34 826 76369.3725 76056.9583 0.41% 908 299s 2820 1434 76287.2839 34 737 76369.3725 76056.9583 0.41% 913 301s 2860 1416 cutoff 35 76369.3725 76062.4630 0.40% 917 306s 2931 1403 76246.5849 32 596 76369.3725 76062.4630 0.40% 910 311s 3102 1425 76363.4387 42 309 76369.3725 76062.4630 0.40% 884 316s 3288 1381 76158.6131 32 741 76369.3725 76063.6886 0.40% 849 321s 3377 1363 76302.7944 37 437 76369.3725 76063.6886 0.40% 845 325s 3436 1347 76346.7969 39 459 76369.3725 76063.6886 0.40% 838 332s H 3437 1274 76369.372471 76063.6886 0.40% 838 332s 3458 1290 76350.3118 40 499 76369.3725 76063.6886 0.40% 835 335s 3694 1243 76342.0630 61 481 76369.3725 76272.0081 0.13% 799 340s 3962 1172 cutoff 58 76369.3725 76336.7477 0.04% 752 345s 4260 1082 76354.0817 71 523 76369.3725 76342.6575 0.03% 707 350s 4666 947 76350.0466 63 589 76369.3725 76346.2908 0.03% 654 355s 5087 802 76355.9659 57 457 76369.3725 76350.6848 0.02% 609 363s 5348 730 76358.2772 60 522 76369.3725 76353.1382 0.02% 582 3685848 496 cutoff 72 76369.3725 76356.8278 0.02% 537 374s 6911 17 76364.4056 72 252 76369.3725 76358.3448 0.01% 462 379s Cutting planes: Gomory: 2 Cover: 2 Projected implied bound: 2 MIR: 101 Mixing: 2 StrongCG: 5 Flow cover: 22 GUB cover: 1 Inf proof: 3 Zero half: 5 RLT: 6 Relax-and-lift: 2 Explored 7738 nodes (3273161 simplex iterations) in 379.66 seconds (519.24 work units) Thread count was 12 (of 12 available processors) Solution count 10: 76369.4 76369.4 76370.7 ... 95199.7 Optimal solution found (tolerance 1.00e-04) Best objective 7.636937247086e+04, best bound 7.636937247086e+04, gap 0.0000%
Please see the .LP file of this constraint and have a look at Spinningusage[125,85] constraint only.
File : https://drive.google.com/file/d/1NkGdJF2BnobqRWsaR4YjC6AAkRANUZz6/view?usp=sharing
Now please take a look at this solution file and see Tasks[17,80]. Its value is 1. Have a look.
File : https://drive.google.com/file/d/1eFpmm8QB-pFL0VW-IlBEOs1m6Yjcvgun/view?usp=sharing
Now look at the SpinningResource[125,85] in this solution file. It should be 34. How come it is 0 ?
It is a simple calculation from these two files.
I am trying to understand this unusual error but unable to find the glitch.
Thank you,
Margi
0 -
-
-
Now look at the SpinningResource[125,85] in this solution file. It should be 34. How come it is 0 ?
Everything is correct. Using the solution values, constraint SpinningResource[125,85] reads
- 34 Tasks[17,80] + SpinningResource[125,85] <= 0
=> - 34 * 1 + 0 <= 0
<=> 0 <= 34which is a correct inequality. So it looks like setting \(\texttt{Tasks[17,80]=1}\) is required somewhere else and \(\texttt{SpinningResource[125,85]}\) is forced to \(0\) from somewhere else.
It is probably the case that this set of constraints
R16926: SpinningResource[125,85] - SpinningResource[125,86] = 0
R16927: SpinningResource[125,85] - SpinningResource[125,87] = 0
R16928: SpinningResource[125,85] - SpinningResource[125,88] = 0
R16929: SpinningResource[125,86] - SpinningResource[125,87] = 0
R16930: SpinningResource[125,86] - SpinningResource[125,88] = 0
R16931: SpinningResource[125,87] - SpinningResource[125,88] = 0
R16932: SpinningResource[125,89] - SpinningResource[125,90] = 0
R16933: SpinningResource[125,89] - SpinningResource[125,91] = 0
R16934: SpinningResource[125,89] - SpinningResource[125,92] = 0
R16935: SpinningResource[125,90] - SpinningResource[125,91] = 0
R16936: SpinningResource[125,90] - SpinningResource[125,92] = 0
R16937: SpinningResource[125,91] - SpinningResource[125,92] = 0
R16938: SpinningResource[125,93] - SpinningResource[125,94] = 0
R16939: SpinningResource[125,93] - SpinningResource[125,95] = 0
R16940: SpinningResource[125,93] - SpinningResource[125,96] = 0
R16941: SpinningResource[125,94] - SpinningResource[125,95] = 0
R16942: SpinningResource[125,94] - SpinningResource[125,96] = 0
R16943: SpinningResource[125,95] - SpinningResource[125,96] = 0forces it to \(0\), despite the fact that \(\texttt{SpinningResource[125,85]}\) would be allowed to be \(>0\).
Best regards,
Jaromił0 -
-
Hi Jaromił Najman,
When I run the MILP, sometimes it shows like below :
See the warning.
Does it mean anything with the optimal result ?56332 15177 100409.197 75 753 100424.350 100348.931 0.08% 514 1957s
H56340 15177 100424.35033 100348.931 0.08% 514 1957s
56690 15117 cutoff 72 100424.350 100349.876 0.07% 515 2006s
H57125 14921 100423.61460 100350.094 0.07% 516 2006s
58185 14781 100398.161 76 830 100423.615 100352.565 0.07% 518 2065s
59681 14574 100411.355 74 1025 100423.615 100356.098 0.07% 522 2112s
61063 14312 100420.766 73 733 100423.615 100359.914 0.06% 529 2171s
62310 14082 cutoff 77 100423.615 100363.616 0.06% 535 2224s
63770 13897 cutoff 71 100423.615 100367.853 0.06% 541 2281s
65040 13490 100414.794 84 850 100423.615 100372.046 0.05% 546 2341s
*65209 13490 98 100423.61446 100372.127 0.05% 547 2341s
66538 13688 cutoff 74 100423.614 100377.932 0.05% 551 2397s
69564 13836 100414.208 89 664 100423.614 100384.277 0.04% 545 2453s
72215 13864 100392.663 75 1010 100423.614 100387.590 0.04% 537 2480s
73878 13860 100405.075 84 697 100423.614 100389.628 0.03% 533 2505s
75356 13863 cutoff 78 100423.614 100391.065 0.03% 530 2532s
76956 13852 100415.380 81 842 100423.614 100392.544 0.03% 527 2567s
78557 13832 100418.006 87 958 100423.614 100394.024 0.03% 523 2600s
79876 13787 cutoff 86 100423.614 100394.920 0.03% 521 2626s
81408 13645 cutoff 76 100423.614 100396.164 0.03% 519 2656s
82806 13512 cutoff 84 100423.614 100397.136 0.03% 517 2692s
84156 13309 cutoff 79 100423.614 100398.134 0.03% 516 2726s
85527 13076 100409.509 77 836 100423.614 100399.188 0.02% 515 2755s
86661 12740 cutoff 83 100423.614 100399.939 0.02% 515 2787s
87832 12369 cutoff 77 100423.614 100401.202 0.02% 516 2828s
88989 11980 cutoff 83 100423.614 100402.298 0.02% 517 2862s
90135 11488 100419.636 98 649 100423.614 100403.313 0.02% 518 2894s
91404 10913 cutoff 74 100423.614 100404.739 0.02% 518 2937s
92817 10213 cutoff 123 100423.614 100405.925 0.02% 519 2985s
94217 9543 cutoff 101 100423.614 100407.958 0.02% 519 3026s
96345 8877 infeasible 88 100423.614 100409.917 0.01% 516 3081s
98856 8055 cutoff 78 100423.614 100411.685 0.01% 511 3132s
101814 7332 100420.024 142 382 100423.614 100413.224 0.01% 504 3171s
Cutting planes:
Gomory: 44
Cover: 231
Implied bound: 1
Projected implied bound: 18
Clique: 2
MIR: 136
Mixing: 2
StrongCG: 8
Flow cover: 77
Inf proof: 8
Zero half: 46
Mod-K: 1
RLT: 19
Relax-and-lift: 63
Explored 105062 nodes (52212531 simplex iterations) in 3172.34 seconds (3034.75 work units)
Thread count was 12 (of 12 available processors)
Solution count 10: 100424 100424 100424 ... 100432
Optimal solution found (tolerance 1.00e-04)
Warning: max constraint violation (1.6190e-06) exceeds tolerance
Warning: max general constraint violation (1.6190e-06) exceeds tolerance
Best objective 1.004236144608e+05, best bound 1.004146711839e+05, gap 0.0089%0 -
-
-
Hi Margi,
The warning means that there is a constraint which has a violation slightly above the default FeasibilityTol parameter value. It is usually not an issue as long as you check the optimal solution point before using it for your application.
In this particular case, some general constraint is violated. Which general constraints do you use in your model? General constraints are the ones which are added via an addGenConstrXxx() method.
You can usually decrease the violation by decreasing the value of the FeasibilityTol parameter or for general constraints, by either tightening the bounds of the argument variables and/or experimenting the parameters PreSOS1BigM and PreSOS2BigM.
Best regards,
Jaromił0 -
-
Hi Jaromił Najman,
Thank you for letting me know about this.
I was just wondering if optimization programs can use the GPU?
def buildModel():
############################################ price and generation data #############################################
####################### Electricity price ##########################################################
energy_price_hour = [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]
# making the list as per 96 time slots
price_energy = []
for price in energy_price_hour:
price_energy = price_energy + [int(price)] * int(60 / rtn_t0)..........................................
;................................Objective_value = []
for i in range(no_of_episodes):
print(f"Running episode {i}")
objval = buildModel().cuda()
Objective_value.append(objval)
print("Objective value:", Objective_value)I try the function to use cuda if possible, but it seems throwing this error.
Is there a way to enable cuda ? I am sure my system has inbuilt GPU, as seen below:
Best regards,
Margi0 -
-
-
Hi Margi,
Currently, it is not possible to execute Gurobi on a GPU.
Best regards,
Jaromił0 -
Please sign in to leave a comment.
Comments
51 comments