• Gurobi Staff

Hi Margi,

You can formulate $$\max$$ terms using the addGenConstrMax method. Note that this method only accepts single optimization variables as inputs. Thus, you will have to introduce auxiliary variables and constraints. In your case you will have to introduce something similar to

\begin{align*} \text{aux}_{res,PV,t} &= \text{EnergyResource}_{res,t} - E_{PV,t} \quad \forall res,PV,t\\ \text{maxaux}_{res,PV,t} &= \max \{0, \text{aux}_{res,PV,t}\} \quad \forall res,PV,t \end{align*}

You can then use variables $$\text{maxaux}_{res,PV,t}$$ in your objective function.

Your code might look something like

aux = m.addVars(reslist2, PVlist, T, lb=-GRB.INFINITY,name="aux")maxaux = m.addVars(reslist2, PVlist, T, name="maxaux")for res in reslist2:  for PV in PVlist:    for t in T:      # add constr aux = EnergyResouce - E      m.addConstr(aux[res,PV,t] == EnergyResource[res,t] - E[PV,t],name="aux_constr_%s_%s_%d"%(res,PV,t))      # add constr maxaux = max(0, aux)      m.addGenConstrMax(maxaux[res,PV,t],[aux[res,PV,t]],0.0,name="maxconstr_%s_%s_%d"%(res,PV,t))

Best regards,
Jaromił

Many thanks for the reply. I understand why you add these auxiliary variables.

However, I have few questions,

1. E_PV=[866.23,789.27,765.65,740.63,686.67,708.26,720.81,690.84,725.12,740.48,765.95,841.81,      898.6,935.81,924.69,913.97,915.17,896.59,863.81,877.96,861.37,878.02,930.13,1046.18,      1158.19,1107.41,1035.81,992.84,1027.11,1072.46,1112.38,1116.97,1184.55,1234.48,1372.79,      1562.12,1596.25,1620.23,1676.42,1677.04,1668.11,1615.57,1418.05,1254.78,1115.44,1042.66,      965.08,925.74,906.32,856.2,814.01,843.43,896.17,927.68,966.43,944,961.58,1027.75,980.66,      952.12,958.43,1000.16,986.89,988.2,923.57,916.67,886.77,935.58,891.83,889.45,891.55,909.11,      957.5,986.71,1005.58,965.28,924.97,946.16,947.49,942.25,885.71,869.43,842.47,828.72,798.17,      847.9,838.99,747.41,725.22,758.07,758.82,760.96,768.38,719.37,693.03,680.75]
2. aux = m.addVars(reslist2, E_PV, time, lb=-GRB.INFINITY,name="aux")maxaux = m.addVars(reslist2, E_PV, time, name="maxaux")
3.     for res in reslist2:    for PV in E_PV:    for t in time:        m.addConstr(aux[res,PV,t] == EnergyResource[res,t] - E[PV,t],name="aux_constr_%s_%s_%d"%(res,PV,t))        m.addGenConstrMax(maxaux[res,PV,t],[aux[res,PV,t]],0.0,name="maxconstr_%s_%s_%d"%(res,PV,t))

My questions are:

1. What is E[PV,t] ?, Should not it be aux[PV,t]

2. Also these are the constraints we are adding, but how to modify my objective function?

3. This equation is bit confusing.

m.addGenConstrMax(maxaux[res,PV,t],[aux[res,PV,t]],0.0,name="maxconstr_%s_%s_%d"%(res,PV,t))

Also, the objective function, Current one:

obj=gp.quicksum(price_energy[t-1]*EnergyResource[res,t] for res in reslist2 for t in time1)

How to modify it?

I am trying to understand these methods, but little confused with the abbreviations you have used, so if you can just detail those?

-Best Regards,

Margi

• Gurobi Staff

Hi Margi,

1. What is E[PV,t] ?, Should not it be aux[PV,t]

$$\texttt{E[PV,t]}$$ refers to your first comment where you show the subtraction $$\texttt{EnergyResource[res,t]- E[PV,t]}$$ in a screenshot. Since in your latest comment $$\texttt{E_PV}$$ is just defined over $$\texttt{t}$$ (I guess?) you have to introduce the auxiliary variables only for res and t

aux = m.addVars(reslist2, time1, lb=-GRB.INFINITY,name="aux")maxaux = m.addVars(reslist2, time1, name="maxaux")

Also these are the constraints we are adding, but how to modify my objective function?

Instead of EnergyResource, you can now use the $$\texttt{maxaux}$$ variables

obj=gp.quicksum(price_energy[t-1]*maxaux[res,t] for res in reslist2 for t in time1)

This equation is bit confusing.

Why exactly is this confusing? The equation

m.addGenConstrMax(maxaux[res,t],[aux[res,t]],0.0,name="maxconstr_%s_%d"%(res,t))

defines the constraint $$\text{maxaux}_{res,t} = \max \{0, \text{aux}_{res,t} \}$$. You have to define it like this due to the required function arguments of the addGenConstrMax method which accepts only single optimization variables. This is also the reason why you have to introduce the auxiliary equality constraints

m.addConstr(aux[res,t] == EnergyResource[res,t] - E_PV[t],name="aux_constr_%s_%d"%(res,t))

which model $$\text{aux}_{res,t} = \text{EnergyResource}_{res,t} - EPV_t$$. Please note that you now also have to adjust the $$\texttt{for}$$-loops to

for res in reslist2:  for t in time1:

Best regards,
Jaromił

Hi Jaromil,

I made the necessary changes as below:

#Adding of auxillary variables for modifying objective function,incuding E_PVaux = steel.addVars(reslist2, time1, lb=-GRB.INFINITY,name="aux")maxaux = steel.addVars(reslist2, time1, name="maxaux")

for res in reslist2:    for t in time1:         # add constr aux = EnergyResouce - E        steel.addConstr(aux[res,t] == EnergyResource[res,t] - E_PV[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))

Modification of objective function:

obj=gp.quicksum(price_energy[t-1]*maxaux[res,t] for res in reslist2 for t in time1)

After Compiling, I got that that the model is infeasible, then I did  steel.computeIIS() and steel.write("steelIIS.ilp"),

\ Model steel RTN_copy\ LP format - for model browsing. Use MPS format to capture full model detail.Minimize Subject To aux_constr_174_42: - EnergyResource[174,42] + aux[174,42] = -1418.05 aux_constr_174_42: - EnergyResource[174,42] + aux[174,42] = -1615.57Bounds EnergyResource[174,42] free aux[174,42] freeEnd

Now I am unsure why aux_constr_174_42 is repeating twice? And what this error actually means? Because we want the max(0,aux[res,t]) , So what's the issue in this formulation?

-Best Regards,

Margi

• Gurobi Staff

Hi Margi,

Now I am unsure why aux_constr_174_42 is repeating twice?

How are $$\texttt{reslist2}$$ and $$\texttt{time1}$$ defined? Do you have duplicate entries in $$\texttt{reslist2}$$? This might cause trouble. You have to find out why the combination $$\texttt{res=174,t=42}$$ occurs twice in the $$\texttt{for}$$-loop to fix the issue.

And what this error actually means? Because we want the max(0,aux[res,t]) , So what's the issue in this formulation?

The model is infeasible because you state that $$x+y = -1418.05$$ and simultaneously $$x + y = -1615.57$$ which is not possible.

Best regards,
Jaromił

Hi Jaromil,

Thanks for pointing it out.I will look into that.

I have one simple question:

aux = steel.addVars(reslist2, time1, lb=-GRB.INFINITY,name="aux")

Why did you put lb=-GRB.INFINITY here?

• Gurobi Staff

Hi Margi,

Why did you put lb=-GRB.INFINITY here?

The default lower bound for variables in Gurobi is $$0$$. But the value of EnergyResource - EPV might be negative. Thus, I set the lower bound of these auxiliary variables to be $$-\infty$$ to allow for negative values.

This is not needed for the max auxiliary variables because you are taking $$\max\{0, \text{aux}\}$$ which is always $$\geq 0$$.

Best regards,
Jaromił

Hi Jaromil,

I see ! Thank you for this wonderful information .

I will continue rectifying the error .

Kind regards
-Margi

Thanks for the information. I managed to somehow rectify the error.

I have one question: How to add custom heuristics in my existing problem for faster MIP's? I came to this video: Faster MIPs Using Custom Heuristics - Gurobi Optimization , but literally I am unaware of how to do so in my existing problem.

I am aware of Heuristics - Gurobi Optimization but I want to check on this : Faster MIPs Using Custom Heuristics - Gurobi Optimization , just to check the problem's performance in presence of heuristics?

What kind of modifications do I need to make in the existing mathematical formulation?

-Kind regards,

Margi

• Gurobi Staff

Hi Margi,

How to add custom heuristics in my existing problem for faster MIP's?

There are multiple ways to adjust heuristics and provide feasible points to Gurobi.

You could compute a feasible point for your model before starting the optimization. You can provide the feasible point to Gurobi through a MIPStart.

You can control Gurobi's built-in heuristics, via the parameters Heuristics, RINS,  and NoRelHeurTime. The Heuristics parameter is a meta parameter which controls the overall effort spent by Gurobi for heuristics. The RINS parameter controls the frequency of the RINS heuristic. The NoRelHeurTime parameter controls the time spent in the "No Relaxation Heuristic". It is a special heuristic which is not run per default. This heuristic is executed before the first root node relaxation is solved and searches for a feasible point without solving the root relaxation (thus the name).

A way to provide feasible points to Gurobi computed by your own heuristics is to use callbacks. You have seen the advanced webinar on this topic already Faster MIPs Using Custom Heuristics - Gurobi Optimization. You might also be interested in the video on Python III: Optimization and Heuristics.

but literally I am unaware of how to do so in my existing problem.

Are you unaware how to compute a feasible point to your model? This is not something that can be answered easily. You would have to use specific problem knowledge and think of a new heuristic tailored to your particular problem to compute a feasible point and provide it to Gurobi through a callback.

Or don't you know how to provide a feasible point to Gurobi? In this case, please refer to How do I use MIP starts? and the videos you have watched already. In particular, the videos on introduction to callbacks.

What kind of modifications do I need to make in the existing mathematical formulation?

You don't have to do any modifications to your existing mathematical formulation, because feasible points to not alter the given formulation.

Best regards,
Jaromił

Hi Jaromil,

Many thanks for this detailed information. I looked upon the video lecture and related documentation from Gurobi. It was bit clear now.

Kind regards,

Margi

Can I ask one question as I am seeing weird error in gurobi results:

I have defined Tasks variables as BINARY:

In results I see :

How it is possible? It should either be 0 or 1.

Can you please let me know?

• Gurobi Staff

This is answered in the Knowledge Base article Why does Gurobi sometimes return non-integral values for integer variables?

Thanks Jaromil.

If I have correctly understood, there is nothing wrong with the optimisation. Its just the rounding feature. Am I correct?

Sorry but just to double check.

• Gurobi Staff

Am I correct?

Correct.