• Gurobi Staff

Hi Johannes,

To be concise, we can use a single binary variable $$z$$ to represent the sign of the variable $$m$$:

\begin{align*}z = 1 &\iff m \geq 0.\end{align*}

We can equivalently write this using the following two implications:

\begin{align*}z = 1 &\implies m \geq 0 \\ z = 0 &\implies m < 0.\end{align*}

Unfortunately, we can't directly model strict inequalities like $$m < 0$$ in a mathematical programming framework. However, we can approximate this inequality with $$m \leq -\epsilon$$ for a sufficiently small $$\epsilon > 0$$ (maybe $$\epsilon = 10^{-6}$$):

eps = 1e-6model.addConstr((z == 1) >> (m >= 0))model.addConstr((z == 0) >> (m <= -eps))

Then, you can change your current indicator constraints to use this single binary variable $$z$$:

model.addConstr((z == 1) >> (t1 - t1a == 0))model.addConstr((z == 0) >> (t2 - t1b == 0))

Thanks,

Eli

Dear Eli,

Thank you very much for your instant response. I implemented your suggested approach in the model. However, it is still not feasible and gurobi says unbounded. I also added some lb and ub where known and set the objective function to zero like suggested in another post. The amount of variables is equal the amount of functions.

I would really like to share the mathematical model with somebody who can help me in finding a solution. It is all about district heating simulation (massflows, temperatures, pressures and heat losses).

Can you further assist me in this topic or do you know somebody who could?

Best regards,

Johannes

• Gurobi Staff

Hi Johannes,

The model must be infeasible due to some other constraints, then. Did you declare a lower bound of $$\texttt{-GRB.INFINITY}$$ when you created the $$\texttt{m}$$ variable? Note that the lower bound of a variable is $$0$$ by default.

To determine the source of the infeasibility, you could try computing an Irreducible Inconsistent Subsystem with Model.computeIIS(). After doing so, write out an .ilp file and visually inspect it to identify which constraints are contributing to the infeasibility. E.g.:

model.computeIIS()model.write('model.ilp')

Feel free to post the contents of the .ilp file here if it's not clear which constraints are causing the infeasibility.

Thanks,

Eli

Dear Eli,

The method model.computeIIS() is not giving any results. Moreover, it is not starting properly. It stucks when calling without error message and the CPU keeps running at full load (test 10 min).

Here is a shortened version of my model. I kicked out all duplicates and tried to keep it to a minimum information. That means the equations are given at least once and also the variables. Of course, any equation and variable is more frequently in the real model. Hope this will help solving the problem. Before gurobi I have tried ipopt (interior point optimization) and the model worked with a good starting value but had trouble to optimize the objfnc.

\ Model sim4dhs
\ LP format - for model browsing. Use MPS format to capture full model detail.
Minimize
0 t[0] + 0 t[1] + 0 t[2] + 0 t[4] + 0 z[0] + 0 z[1] + 0 z[2] + 0 z[3]
+ 0 z[4] + 0 z[5]
Subject To
mbs[0]: - m_pipes[0] + m_cons[2] = 0

eb2_cons[0]: - t[5] + ta[6] = 0
eb2_prod[0]: - t[3] + ta[9] = 0

ib1[0]: - p[0] + pa[0] = 0

ib2[0]: p[1] - pb[0] = 0

ARG1Q[0]: - q_pipes[0] - ARG1[0] + ARG1Q[0] = 0

ppd_drop_lin[0]: - 2431.708407416106 m_pipes[0] + pa[0] - pb[0] = -29430

cons_q[0]: - q_cons[0] = 9999.999999971065

cons_tb[0]: - tb[6] = -313.15
prod_tb[0]: - tb[9] = -373.15

eb1[0]: [ - m_pipes[0] * ta[0] + m_cons[2] * tb[8] ] = 0

mfs[0]: q_pipes[0] + [ 4186.8 m_pipes[0] * ta[0]
- 4186.8 m_pipes[0] * tb[0] ] = 0

mfs[6]: q_cons[0] + [ 4186.8 m_cons[0] * ta[6] - 4186.8 m_cons[0] * tb[6]
] = 0
mfs[9]: q_prod[0] + [ 4186.8 m_prod[0] * ta[9] - 4186.8 m_prod[0] * tb[9]
] = 0

cnstrARGa10: ARG1[0] + 1.1771188200000001e+06 absM[0] + [
- 4186.8 tab[0] * absM[0] ] = 0

pipes_q[0]: [ - logARG1[0] * absM[0] + logARG1Q[0] * absM[0] ]
= -0.0450213908511005

flowdirectpos[0]: z[0] = 1 -> m_pipes[0] >= -0

flowdirectneg[0]: z[0] = 0 -> m_pipes[0] <= -1e-06

eb2_pipes_pos[0]: z[0] = 1 -> - t[0] - 0 t[1] - 0 t[2] - 0 t[3] - 0 t[4]
- 0 t[5] - 0 t[6] - 0 t[7] + ta[0] = -0

eb2_pipes_neg[0]: z[0] = 0 -> - 0 t[0] - t[1] - 0 t[2] - 0 t[3] - 0 t[4]
- 0 t[5] - 0 t[6] - 0 t[7] + tb[0] = -0

GC24: z[0] = 1 -> tab[0] - ta[0] = -0

GC30: z[0] = 0 -> tab[0] - tb[0] = -0

Bounds
m_pipes[0] free

-infinity <= q_pipes[0] <= 0

-9999.999999971065 <= q_cons[0] <= 0

281.15 <= t[0] <= 373.15

281.15 <= ta[0] <= 373.15

281.15 <= tb[0] <= 373.15

281.15 <= tab[0] <= 373.15

Binaries
z[0] z[1] z[2] z[3] z[4] z[5]
General Constraints
logARG10: logARG1[0] = LOG ( ARG1[0] )

logARG1Q0: logARG1Q[0] = LOG ( ARG1Q[0] )

absM0: absM[0] = ABS ( m_pipes[0] )

End

model.optimize() starts properly but does not find a solution:

Gurobi Optimizer version 9.0.2 build v9.0.2rc0 (win64)
Optimize a model with 50 rows, 118 columns and 109 nonzeros
Model fingerprint: 0x771caae4
Model has 54 general constraints
Variable types: 112 continuous, 6 integer (6 binary)
Coefficient statistics:
Matrix range [1e+00, 2e+03]
QMatrix range [1e+00, 4e+03]
QLMatrix range [1e+00, 1e+06]
Objective range [0e+00, 0e+00]
Bounds range [1e+00, 4e+04]
RHS range [3e+02, 4e+04]
QRHS range [5e-02, 5e-02]
Warning for adding constraints: zero or small (< 1e-13) coefficients, ignored
Presolve added 0 rows and 2146 columns
Presolve removed 12 rows and 0 columns
Presolve time: 0.01s
Presolved: 216 rows, 2293 columns, 7078 nonzeros
Presolved model has 12 SOS constraint(s)
Presolved model has 40 bilinear constraint(s)
Variable types: 2293 continuous, 0 integer (0 binary)
Root relaxation: objective 0.000000e+00, 75 iterations, 0.00 seconds
Nodes | Current Node | Objective Bounds | Work
Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time
0 0 0.00000 0 40 - 0.00000 - - 0s
0 0 0.00000 0 34 - 0.00000 - - 0s
0 0 0.00000 0 34 - 0.00000 - - 0s
0 0 0.00000 0 34 - 0.00000 - - 0s
0 0 0.00000 0 34 - 0.00000 - - 0s
0 0 0.00000 0 34 - 0.00000 - - 0s
0 2 0.00000 0 45 - 0.00000 - - 0s
9385 2286 infeasible 70 - 0.00000 - 8.0 5s
23343 3650 0.00000 73 37 - 0.00000 - 6.2 10s
36229 4061 infeasible 49 - 0.00000 - 5.7 15s

...and so forth...

Best regards,

Johannes

• Gurobi Staff

Hi Johannes,

I'm a bit confused. Does Gurobi ever declare the model infeasible? The log snippet you posted shows that Gurobi has not found a solution yet, but that does not necessarily mean that the model is infeasible. Note that calculating an IIS can be a challenging problem in itself, though there isn't much reason to try this unless Gurobi terminated with an infeasible status code.

It looks like you have a difficult problem. Problems with nonconvex quadratic constraints are generally hard to solve, even without considering additional integer variables and general constraints. To look for a solution faster, you could try setting MIPFocus=1 or increasing the value of the Heuristics parameter. You could also try passing an initial solution (that you or Ipopt calculate) to Gurobi by setting the Start attribute for each variable.

Thanks,

Eli