Indeed, I was computing the runtime of the dual problem in Benders decomposition. I defined the dual model outside the callback function. In each iteration, I redefined the dual's objective function and then optimized it. I printed the dual's Runtime every time. But it was strange that the runtime could be up to about 20s many times.

• Gurobi Staff

Could you please provide some log output or even better, a minimal working example reproducing this issue?

It might be possible that a previous solution which is used in the second optimization run has a negative impact on the optimization path of the second run and it is indeed better to just re-run the second model from scratch.

Best regards,
Jaromił

I define the sub_model's variables and constraints outside the callback function

In my callback function:

### Retrieve the current value of the master problemr_1 = np.array([[model.cbGetSolution(model.getVarByName('r_' + str(i) + '_' + str(j))) for j in range(N + 1)] for i in range(N + 1)])z_1 = np.array([[[model.cbGetSolution(model.getVarByName('z_' + str(i) + '_' + str(j) + '_' + str(l))) for l in range(L)] for j in range(N + 1)]for i in range(N + 1)])### Retrieve the current value of master problem as the lower boundLB = model.cbGet(GRB.Callback.MIPSOL_OBJBND)
### the w and u variable are the dual variables while the z_1 is the given primal variable fed by the outer problemm_sub.setObjective(    quicksum(w_3[i, j] * y[i - 1, j] for j in range(len(B)) for i in np.arange(N) + 1)    - quicksum(w_4[i, j] * y[i - 1, j] for j in range(len(B)) for i in np.arange(N) + 1)    + quicksum(w_6[i] * D[i - 1, 0] for i in np.arange(N) + 1)    - quicksum(w_7[i] * D[i - 1, 1] for i in np.arange(N) + 1)    - K    * quicksum(        z_1[i][j][l] * u[i, j, l]        for l in range(L)        for j in range(N + 1)        for i in range(N + 1)    )    - -K * quicksum(w_8[i] for i in np.arange(N) + 1),    GRB.MAXIMIZE,)### Print the sloving time at each callback usagestart = time.time()m_sub.optimize()end = time.time()print(m_sub.Runtime)print(end-start)

Below is some part of the result:

Root relaxation: objective 4.700000e+01, 155 iterations, 0.00 seconds (0.00 work units)    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0   47.00000    0   28          -   47.00000      -     -    0s     0     0   47.00000    0   52          -   47.00000      -     -    0s     0     0   47.50000    0   16          -   47.50000      -     -    0s     0     0   47.50000    0   16          -   47.50000      -     -    0s     0     2   48.00000    0   20          -   48.00000      -     -    0s0.158318996429443360.207134962081909186.3162140846252446.322864294052124   105    59   48.00000    7    -          -   48.00000      -   4.9    8s9.4442200660705579.454262018203735   301   155   66.00000   12   34          -   48.00000      -   6.0   18s15.08434987068176315.092468976974487   537   298   70.00000    6   10          -   48.00000      -  11.2   34s9.4082179069519049.416763067245483   910   565   70.00000    9   14          -   48.00000      -  10.6   44s23.0749781131744423.083630084991455

It's abnormal that I should solve the LP dual problem in 10+ seconds. So I included the entire dual model into the callback function(redefine the entire dual model at each time when entering the callback function). All the dual were solved in 1 second.

Since only the z_1 changed in the dual objective function, I changed the objective function's coefficient partly instead of redefining the entire objective function. But it converged to a wrong value. I don't know how to explain it。。。

for i in range(N + 1):    for j in range(N + 1):        for l in range(L):            u[i, j, l].Obj = -K * z_1[i][j][l]
• Gurobi Staff

I am not sure if I fully understand the situation. You are saying that when you solve a separate problem S in a callback of a master model M, then S gets solved faster than when S is solved as a standalone problem outside of M. Is this correct? Could you try to write an MPS/LP file of S and solve it outside of M to rule out any inconsistencies. Depending on how you construct S, it is possible that it uses information available for M, e.g., solution point/basis, which can significantly speed up the solution process.

Best regards,
Jaromił

I think you understood me in an opposite way.

In total, I have 2 ways to define the sub dual model need to be solved in the callback function.

I can define the entire dual model in the callback every time when I enter in callback function.

Or I can define the dual model except the objective function outside the callback function, and only define the objective function in the callback function.

The result shows that for the solve time of a single sub dual variable, it will takes more time if I use the second method.(only change the objective function in the callback function)

The result shows that for the solve time of a single sub dual model, it will takes more time if I use the second method.(only change the objective function in the callback function)

For the total solving process, the first one takes more time. And I found the first method solve more Sub dual problems than the second one( I used a global variable to record its solve times)

• Gurobi Staff

It makes sense that your first method takes more time, because when you rebuild your submodel from scratch every time, you lose solution information, which may speed up the optimization process in a subsequent solve.

I don't know why building the constraints outside of the callback function should make a difference in you solve only a single sub dual model. How significant is the optimization time difference? Could you post two log outputs when solve a single sub dual model? One log output when you construct the model outside of the callback and one output when you construct the complete model inside of the callback function?

Best regards,
Jaromił

Hi, Dr Najman. I found that the sub_dual model didn't presolve when I only redefine its objective function in the callback function. That was why the sub_dual problem would take so much time. Outside the callback function, I set the sub_dual problem's presolve parameter as -1. But the log still showed that it didn't preslove in the callback function.

How should I fix it?

Thanks a lot.

• Gurobi Staff

Did you try explicitly setting the Presolve parameter inside your callback before solving the sub_dual model? You could also try calling the reset method before optimizing and setting the parameter.

Could your provide a minimal reproducible example or at least some code snippet of you are are currently solving the sub_dual model within the callback?

def callBackFunction(model, where):    if where == GRB.Callback.MIPSOL:        ### retrieve the current integer value        r_1 = np.array([[model.cbGetSolution(model.getVarByName('r_' + str(i) + '_' + str(j))) for j in range(N  + 1)] for I in range(N + 1)])        z_1 = np.array([[[model.cbGetSolution(model.getVarByName('z_' + str(i) + '_' + str(j) + '_' + str(l))) for l in range(L)] for j in range(N + 1)] for i in range(N + 1)])        ### Retrieve the current value of master problem as the lower bound        LB = model.cbGet(GRB.Callback.MIPSOL_OBJBND)        ### every time I only redefine the dual problem's objective function in callback function        m_sub.setObjective(quicksum(w_3[i, j] * y[i - 1, j] for j in range(len(B)) for i in np.arange(N) + 1) - quicksum(w_4[i, j] * y[i - 1, j] for j in range(len(B)) for i in np.arange(N) + 1) + quicksum(w_6[i] * D[i - 1, 0] for i in np.arange(N) + 1) - quicksum(w_7[i] * D[i - 1, 1] for i in np.arange(N) + 1) + - K * quicksum(z_1[i][j][l] * u[i, j, l] forlinrange(L) forjinrange(N + 1) foriinrange(N + 1)) - K * quicksum(w_8[i] for i in np.arange(N) + 1) , GRB.MAXIMIZE)                m_sub.Params.Presolve = -1        start = time.time()        m_sub.optimize()        end = time.time()        print(m_sub.Runtime)        print(end-start)

Even though I reset the presolve parameter in the callback function(Method 2), it would take 10+ S to solve the LP problem (except the first time to solve the dual problem).

It is strange because:

(1) When I retrieved a dual model during the process and solved it through read_data, it only took no more than 1S to solve it.

(2) When I redefined the entire sub dual model in the callback function(Method 1), it only took no more than 1S to solve it.

Another question is:

Method 2 could solve the entire problem more quickly. It showed that the sub dual problem was solved fewer times under method 2. Does it make sense?

• Gurobi Staff

Thank you for posting a snippet of your callback and the clarifications. I think we are slowly getting there.

What might be causing the slowdown is that Gurobi threads interfere with each other. Could you please try limit the Threads parameter for both the master problem and the dual problem solved in each callback. For example, if you have 4 threads, then try assigning 2 threads to the master problem and 2 to the dual problem.

Method 2 could solve the entire problem more quickly. It showed that the sub dual problem was solved fewer times under method 2. Does it make sense?

Yes, it does make sense. It is not unusual that the barrier algorithm (Method=2) solves a model faster than the primal (Method=0) or dual (Method=1) simplex algorithm. If you observe that setting Method=2 benefits the performance of your models, then just set it explicitly.

Did you consider avoiding callbacks? You could solve the master problem with SolutionLimit=1, which would terminate the optimization process as soon as a feasible solution is found. You can then access all necessary information via model attributes, adjust the dual model and solve it. Next, you would add the newly computed constraint to the master problem and resolve it. This whole process would have to be looped.

Okay. Thank you for your quick response!!!! I am trying it(:

Hi, Jaromil! I have tried your suggestion.  However, there are still some issues I want to talk about with you.

m_sub.Params.Threads = 4m_master.Params.Threads = 4
I have 8 threads in total. After assigning, the solving time didn't change too much.

2. "Method"
Maybe I have made some confusion.  It doesn't mean the "Method" parameter.

Way 1: redefine the entire sub-dual model in the callback function
Way 2: only redefine the objective function in the callback function

Since the subdual model contains 106240 rows, 15168 columns, and 313760 nonzeros, it will take much time to redefine them every time.

The result shows that Way 2 is more time-efficient. The sub-dual problem was solved fewer times under Way 2. Thus, though it took more time to solve LP, the total solving time is less. Does it make sense about the fact that the sub-dual problem was solved fewer times under Way 2?

3. Other ways

The model is a complicated MIP model. In the callback function, I deal with the sub-tour elimination in the first part and add benders cut in the second part. All cuts are added in a lazy-constraint way.  I don't know whether the loop is time-efficient...

Thank you so much(:

• Gurobi Staff

Maybe I have made some confusion.  It doesn't mean the "Method" parameter.

Way 1: redefine the entire sub-dual model in the callback function
Way 2: only redefine the objective function in the callback function

Thank you for clarifying. It makes sense now.

Did you try changing objective coefficients via the variable Obj attribute instead of reconstructing the objective from scratch?

Does it make sense about the fact that the sub-dual problem was solved fewer times under Way 2?

Yes, if you only change the objective, then Gurobi can warm-start from a previous solution. This is not possible if you reconstruct the model from scratch. It is also possible that through the warm-start Gurobi finds a slightly different optimal solution point (if, e.g., multiple solutions are available) and ultimately, under Way 2 fewer problems have to be solved.

3. Other ways

The model is a complicated MIP model. In the callback function, I deal with the sub-tour elimination in the first part and add benders cut in the second part. All cuts are added in a lazy-constraint way.  I don't know whether the loop is time-efficient...

Hard to tell whether it would be worth trying a priori. You could keep this in mind in case that nothing else helps.

Did you try changing objective coefficients via the variable Obj attribute instead of reconstructing the objective from scratch?

I have tried before. The results are the same. You mentioned "the warm start". Due to it, the sub-dual LP problem should be solved more quickly. But it took much more time to solve compared with that under Way 1...

• Gurobi Staff

Due to it, the sub-dual LP problem should be solved more quickly. But it took much more time to solve compared with that under Way 1...

This is also possible in some cases. Unfortunately, in general a warm-start does not guarantee performance improvement.

Did you try explicitly setting the Presolve parameter inside your callback before solving the sub_dual model? You could also try calling the reset method before optimizing and setting the parameter.

It doesn't work if I set the Presolve parameter inside your callback before solving the sub_dual model. But there is something wrong if I use reset method(=0, =1). If I reset,  it took normal time to solve the LP problem. But it turned out that the result is wrong. Is there any other method to set the parameter(Through environment)?

• Gurobi Staff

But there is something wrong if I use reset method(=0, =1). If I reset,  it took normal time to solve the LP problem. But it turned out that the result is wrong.

What exactly do you mean by wrong result? Could you please share log files of the runs with and without call of the reset method?

Is there any other method to set the parameter(Through environment)?

You can use the Env.setParam() method.

To be more specific, I use the reset method in the callback function

m_sub.reset(1)# m_sub.Params.Presolve = -1start = time.time()m_sub.optimize()end = time.time()print(m_sub.Runtime)print(m_sub.objVal)print(end-start)### Print the iteration numbernumber = number + 1print(number)

I set the parameter as 0 and 1 separately.  [A value of 1 discards additional information that affects the solution process but not the actual model (currently MIP starts, variable hints, branching priorities, lazy flags, and partition information). The default value 0 just discards the solution.]

The result shows that the solving time of dual problem is normal if I reset model. But the iteration times are much larger if I reset the model which leads to the unacceptable solving time in total. I have mentioned that even though the solving time of sub-dual problem was abnormally large but the iteration times were few.

Does it make sense?

• model.reset(0)
# solving time0.054704904556274414# objective value24.199999999999903# solving time
0.06278800964355469# iterative times
10.048840999603271484 24.199999999999875 0.05006289482116699 20.04844188690185547 24.800000000000143 0.04964804649353027 30.049054861068725586 24.199999999999903 0.05057096481323242 40.04824709892272949 24.199999999999903 0.04944610595703125 5...
• model.reset(1)
0.05489611625671387
24.199999999999903
0.06413078308105469
10.0494990348815918 24.199999999999875 0.05093812942504883 20.04975700378417969 24.800000000000143 0.051056861877441406 30.05099916458129883 24.199999999999903 0.05276298522949219 40.04653191566467285 24.199999999999903 0.047812938690185555
• model(not reset)
0.0669231414794921924.1999999999999030.0810890197753906210.516111135482788124.1999999999998750.520132780075073220.766838073730468824.800000000000110.772233009338378930.635078191757202124.1999999999998750.640252828598022541.36330485343933124.1999999999998751.36778807640075685...(Total times: 11)

Does it make sense? I feed the model different data but it always shows that once i reset the model in the callback, the iterative times will much larger and the model couldn't be solved under limited time. But I only abandoned the sub-dual problem's previous solution. Why would that happen?

• Gurobi Staff

For a better analysis, could you please try to extract the LP files via the write function and try to reproduce this behavior outside of your callback. Extracting the first sub model should be enough.

m_sub.reset(1)# m_sub.Params.Presolve = -1start = time.time()m_sub.optimize()end = time.time()m_sub.write("iteration%d.lp"%(number))print(m_sub.Runtime)print(m_sub.objVal)print(end-start)### Print the iteration numbernumber = number + 1print(number)

You can then try to reproduce the behavior outside of the callback via

m_sub = gp.read("iteration0.lp")# extract variables by name which you need laterstart = time.time()m_sub.optimize()end = time.time()print(end-start)# perform model changesm_sub.reset(1) # or not resetstart = time.time()m_sub.optimize()end = time.time()print(end-start)

This would make it easier to grasp what is happening here. In best case, you could also share the sub model lp file as described in Posting to the Community Forum.

I feed the model different data but it always shows that once i reset the model in the callback, the iterative times will much larger and the model couldn't be solved under limited time.

Could you please share the log output of the sub model optimization? The first 2 iterations should be enough. You can write to log files via the LogFile parameter.

m_sub.setParam("LogFile", "logfile_iteration%d.log"%(number))

Does it make sense?

There are cases where re-using previous solution information may deteriorate optimization performance. However, in your case, the downside of not getting rid of previous solution information is very large (a factor of almost 10x). Thus, we should definitely investigate further as mentioned above in this comment.

On a separate note: Could you please try defining a separate environment for your main model and your sub model? Maybe there are some parameters being reset which were set for the master model.

I feed a data which can be solved by reset method. In the first model, the sub dual problem has been solved for 231 times(6.3411 s in total) to get the final result, which is same to the iteration times(60.9841 s in total) if I define the entire sub dual model in the callback every time vs 12 times (1.3453 s in total) if I don't reset the model.

I extracted the iteration 0's file and ran them as your suggestions.

Set parameter UsernameAcademic license - for non-commercial use only - expires 2023-02-01Read LP format model from file iteration0.lpReading time = 0.02 seconds: 7360 rows, 1398 columns, 20840 nonzerosGurobi Optimizer version 9.5.1 build v9.5.1rc2 (mac64[x86])Thread count: 4 physical cores, 8 logical processors, using up to 8 threadsOptimize a model with 7360 rows, 1398 columns and 20840 nonzerosModel fingerprint: 0x052f33afCoefficient statistics:  Matrix range     [1e+00, 1e+00]  Objective range  [1e+00, 3e+01]  Bounds range     [3e+01, 3e+01]  RHS range        [2e-01, 2e-01]Presolve removed 7105 rows and 1124 columnsPresolve time: 0.01sPresolved: 255 rows, 274 columns, 724 nonzerosIteration    Objective       Primal Inf.    Dual Inf.      Time       0    5.4000000e+32   1.560000e+32   5.400000e+02      0s     158    1.2800000e+01   0.000000e+00   0.000000e+00      0sSolved in 158 iterations and 0.02 seconds (0.01 work units)Optimal objective  1.280000000e+010.0238339900970459Discarded solution information including additional informationGurobi Optimizer version 9.5.1 build v9.5.1rc2 (mac64[x86])Thread count: 4 physical cores, 8 logical processors, using up to 8 threadsOptimize a model with 7360 rows, 1398 columns and 20840 nonzerosModel fingerprint: 0x052f33afCoefficient statistics:Matrix range     [1e+00, 1e+00]  Objective range  [1e+00, 3e+01]  Bounds range     [3e+01, 3e+01]  RHS range        [2e-01, 2e-01]Presolve removed 7105 rows and 1124 columnsPresolve time: 0.00sPresolved: 255 rows, 274 columns, 724 nonzerosIteration    Objective       Primal Inf.    Dual Inf.      Time       0    5.4000000e+32   1.560000e+32   5.400000e+02      0s     158    1.2800000e+01   0.000000e+00   0.000000e+00      0sSolved in 158 iterations and 0.01 seconds (0.01 work units)Optimal objective  1.280000000e+010.006536960601806641

However, in your case, the downside of not getting rid of previous solution information is very large (a factor of almost 10x).

The point lies on why getting rid of previous solution would lead to such larger iteration times? Too many iterations make the model intractable.

Thank you so much!

• Gurobi Staff

Could you also please share the exact code you used to get your posted output?

If I understand correctly, you did not perform any change to the model between the 2 optimizations. But in your callback, you change some coefficients. Could you please also incorporate this model change after the first iteration into your snippet to reproduce the issue?

The point lies on why getting rid of previous solution would lead to such larger iteration times?

I would guess that your algorithm uses the solution point found by the sub model. If you get rid of previous solutions, it is very possible that the Simplex algorithm finds a slightly different solution point with the same solution value. This then may lead to many more iterations needed for your algorithm. You could check this thesis by checking for the solution point in both cases and look for differences.

Could you please also incorporate this model change after the first iteration into your snippet to reproduce the issue?

The necessary code(generate the data): https://drive.google.com/file/d/1wj7epny1yRNouQUzET5Q1zgjiTQU_4eE/view?usp=sharing

For the second part, it's reasonable that the disparity in iteration times results from the multiple solutions to the dual problem. But I tried different scales of data, it always shows that one I reset the model,  the iteration times would explode.

• Gurobi Staff

Thank you for the code. There is a file called $$\texttt{c202.txt}$$ missing.

Sorry for that. Here is the c202.txt.

Thank you!

• Gurobi Staff

Hi Jacob,

m_sub.reset(1) # or 0

the values for $$\texttt{r_1}$$ and $$\texttt{z_1}$$ variables are different from when running with the reset method. This explains the different behavior of the overall algorithm and needing a different number of iterations.

You can verify this by printing $$\texttt{r_1}$$ and $$\texttt{z_1}$$ in each iteration and comparing them.

The reason why each sub model solves faster with reset is that presolve reduces the size of the model and makes the overall solution process quicker. With warm start (without reset), the original model is solved (no presolve is performed). You can change this behavior by setting the LPWarmStart parameter to 2 (use warm start for presolved model - this does degrade performance in your case) or 0 (equivalent to reset).

You could try adding additional constraints to your sub model to lead the overall algorithm into the right direction and ultimately try to reduce the number of iterations of the overall algorithm. I am not an expert on Benders decomposition so i cannot really help here but maybe you are able to find something in the literature.

I hope the above explanation helps.

Best regards,
Jaromił