# Gurobi gets stuck after "Loaded MIP start from previous solve"

**URL:** <https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406>\
**Category:** Optimization (Mathematical)\
**Tags:** question, gurobi, mip\
**Created:** [July 2, 2025, 9:11am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406 "2025-07-02T09:11:34Z")\
**Posts on this page:** 14\
**Page:** 1

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 2, 2025, 9:11am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/1 "2025-07-02T09:11:34Z")

</div>

I’m writing the Benders’ decomposition for an MIP and encountering a very weird behavior of Gurobi.

There is an MIP that I need to modify its obj and repeatedly re-solve (Gurobi can (according to its Logging) `Loaded MIP start from previous solve with objective`). Some obj coefficient might leads to a harder model (surely).

But I discover a model that once become very very slow to resolve. Therefore I set a TimeLimit and intercept its status which at that moment is questionable. Then I write it as an MPS file, and opened a new Julia REPL to directly solve that—it suddenly becomes very fast (at the supposed speed).

The code to generate the questionable model is here

```julia
import JuMP.value as ı
import JuMP, Gurobi; GRB_ENV = Gurobi.Env();
import LinearAlgebra.⋅ as ⋅
import Random

seed = 4599192817953253712 # may need to cherry-pick
Random.seed!(seed)

J = 9; # number of resources
W = rand(9:49, J); # the sole width type of each resource
o = rand(30:0.01:70, J); # occupied cost
u = 7 * rand(0.75:0.01:1.25, J)/30 .* W; # unit price
L = rand(6:17, J); # purchasing num limit
I = 10; # number of customer
w = rand(7:0.7:24, I); # width needed
S = 2; # number of scenes
d = rand(6:35, I, S); # number needed

function Min_s_prim_model(s) # evaluate the 2nd stage exactly ⚠️ Feasible Region is NOT fixed
    Min_s = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV))
    JuMP.@variable(Min_s, m2[1:J]) # JuMP.fix.(m2, m) # 🟡 [Jit] Linking 
    JuMP.@variable(Min_s, 0 ≤ n[1:I, 1:J], Int)
    JuMP.@constraint(Min_s, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
    JuMP.@constraint(Min_s, [j = 1:J], m2[j]W[j] ≥ w ⋅ view(n, :, j))
    JuMP.@objective(Min_s, Min, 0/S)
    return Min_s
end; function Rhs_model_s(s) # 🟢 Feasible Region is fixed
    Rhs = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV));
    JuMP.@variable(Rhs, 0 ≤ b[1:J] ≤ 1, Int);
    JuMP.@variable(Rhs, 0 ≤ m[1:J], Int);
    JuMP.@constraint(Rhs, [j = 1:J], L[j]b[j] ≥ m[j]);
    JuMP.@objective(Rhs, Min, 0) # to be filled
    JuMP.@variable(Rhs, 0 ≤ n[1:I, 1:J], Int)
    JuMP.@constraint(Rhs, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
    JuMP.@constraint(Rhs, [j = 1:J], m[j]W[j] ≥ w ⋅ view(n, :, j));
    ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### 
    # Rhs_1 = Rhs_model_vec[1];
    # μ = -0.1 * o
    # π = 0.1 * u
    # JuMP.set_objective_coefficient.(Rhs_1, Rhs_1[:b], μ)
    # JuMP.set_objective_coefficient.(Rhs_1, Rhs_1[:m], π)
    # JuMP.optimize!(Rhs_1); JuMP.assert_is_solved_and_feasible(Rhs_1, allow_local = false)
    # JuMP.objective_value(Rhs_1)
    ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### 
    return Rhs
end; function Rhs_surrogate_plus_s(s) # 🔵 θ-Cuts can be reused, whereas `varycon` needs to be modified
    Sur = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV));
    JuMP.@variable(Sur, ϕ ≤ 500.0); # ⚠️ initial upper bound need to tune
    JuMP.@objective(Sur, Max, ϕ)
    JuMP.@variable(Sur, θ); # This is the `Rhs_surrogate`, And the follows are `plus`
    JuMP.@variable(Sur, μ[1:J]);
    JuMP.@variable(Sur, π[1:J]);
    # θ-Cuts to be added, which reads `θ ≤ b ⋅ μ + m ⋅ π` in which (b, m) are vertices of the 🟢 Feasible Region of Rhs_model_s(s)
    JuMP.@constraint(Sur, varycon, ϕ - θ ≤ false); # TBD
    return Sur
end

prim_model_vec = [Min_s_prim_model(s) for s = 1:S]; # ✅ evaluate the 2nd stage ObjVals
Rhs_model_vec = [Rhs_model_s(s) for s = 1:S]; # ✅ evaluate the RHS of Lag cut
Sur_model_vec = [Rhs_surrogate_plus_s(s) for s = 1:S]; # ✅ used to generate trial coefficients of Lag cut
JuMP.set_silent.(prim_model_vec); JuMP.set_silent.(Rhs_model_vec); JuMP.set_silent.(Sur_model_vec); 
Mst = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # ✅ used to generate trial 1st-stage solutions
JuMP.@variable(Mst, 0 ≤ b[1:J] ≤ 1, Int);
JuMP.@variable(Mst, 0 ≤ m[1:J], Int);
JuMP.@constraint(Mst, [j = 1:J], L[j]b[j] ≥ m[j]);
JuMP.@objective(Mst, Min, o ⋅ b + u ⋅ m);
JuMP.set_silent.(Mst); 

JuMP.unset_silent(Rhs_model_vec[1]) # This might be hard
JuMP.set_attribute(Rhs_model_vec[1], "TimeLimit", 60)
suspect = nothing

JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
B = ı.(b); # to be globally ✂️ off
M = ı.(m); # to be globally ✂️ off
while true # iteratively find a violating cut for the Master problem
    infeasible_block = 0
    for s = 1:S
        Min_s = prim_model_vec[s]; JuMP.fix.(Min_s[:m2], M);  
        JuMP.optimize!(Min_s); termination = JuMP.termination_status(Min_s)
        if termination == JuMP.INFEASIBLE
            infeasible_block = s
            break # [early break]
        elseif termination != JuMP.OPTIMAL
            error("unknown termination status")
        end
    end
    s = infeasible_block;
    s == 0 && error("We have now derived a 1st-stage FEASIBLE solution!")
    Sur = Sur_model_vec[s]; varycon = Sur[:varycon];
    θ, μ, π = Sur[:θ], Sur[:μ], Sur[:π];
    JuMP.set_normalized_coefficient.(varycon, μ, B);
    JuMP.set_normalized_coefficient.(varycon, π, M);
    JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false); ub = JuMP.objective_bound(Sur)
    Μ, Π = ı.(μ), ı.(π); # to be locally ✂️ off
    Θ = ı(θ); # [🪜Upper] to be locally ✂️ off
    rhs = NaN; # [🪜Lower] pre-allocate
    while true # try generate the coefficient of Lag cut until success
        Rhs_s = Rhs_model_vec[s]
        local b, m = Rhs_s[:b], Rhs_s[:m]
        JuMP.set_objective_coefficient.(Rhs_s, b, Μ)
        JuMP.set_objective_coefficient.(Rhs_s, m, Π)
        JuMP.optimize!(Rhs_s); # JuMP.assert_is_solved_and_feasible(Rhs_s, allow_local = false)
        if JuMP.termination_status(Rhs_s) == JuMP.TIME_LIMIT
            suspect = Rhs_s
            return # 🍅🍅🍅
        end
        rhs = JuMP.objective_value(Rhs_s)
        lb = rhs - Μ ⋅ B - Π ⋅ M
        if lb > 4e-3 # 🟣 global COT
            @info "Find a cut for master with vio = $lb"
            break # go back to add cut for the outer Master
        end
        local_B, local_M = ı.(b), ı.(m) # vertices
        vertical_distance = Θ - (local_B ⋅ Μ + local_M ⋅ Π) # 🪜 having this since you approximate `Rhs_model_s`
        if vertical_distance > 8e-4 # 🟣 local COT
            JuMP.@constraint(Sur, θ ≤ local_B ⋅ μ + local_M ⋅ π)
        else
            error("weird 568435686586442")
        end
        JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false); ub = JuMP.objective_bound(Sur)
        Μ, Π = ı.(μ), ı.(π); # to be locally ✂️ off
        Θ = ı(θ); # [🪜Upper] to be locally ✂️ off
    end
    JuMP.@constraint(Mst, Μ ⋅ b + Π ⋅ m ≥ rhs)
    JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
    M = ı.(m); # to be globally ✂️ off
    B = ı.(b); # to be globally ✂️ off
end

Gurobi.GRBwrite(JuMP.unsafe_backend(suspect), "data/suspect.mps")

```

On my machine it branches enormously (If I don’t set a TIME\_LIMIT, it will branch more and non-stop in a reasonable time):

```julia
Gurobi Optimizer version 12.0.2 build v12.0.2rc0 (win64 - Windows 11.0 (26100.2))

CPU model: 11th Gen Intel(R) Core(TM) i5-1135G7 @ 2.40GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads

Non-default parameters:
TimeLimit 60

Optimize a model with 28 rows, 108 columns and 207 nonzeros
Model fingerprint: 0x8bbf69b7
Variable types: 0 continuous, 108 integer (0 binary)
Coefficient statistics:
  Matrix range [1e+00, 5e+01]
  Objective range [2e+01, 3e+04]
  Bounds range [1e+00, 1e+00]
  RHS range [6e+00, 4e+01]

Loaded MIP start from previous solve with objective -1.36424e-12

Presolve removed 7 rows and 7 columns
Presolve time: 0.00s
Presolved: 21 rows, 101 columns, 193 nonzeros
Variable types: 0 continuous, 101 integer (2 binary)

Root relaxation: objective -3.505153e+02, 53 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 -350.51533 0 11 -0.00000 -350.51533 - - 0s
H 0 0 -51.7171719 -350.51533 578% - 0s
     0 0 -339.13976 0 20 -51.71717 -339.13976 556% - 0s
H 0 0 -53.4817681 -338.94230 534% - 0s
     0 0 -327.61043 0 20 -53.48177 -327.61043 513% - 0s
     0 0 -327.61043 0 20 -53.48177 -327.61043 513% - 0s
     0 0 -327.61043 0 19 -53.48177 -327.61043 513% - 0s
     0 0 -327.61043 0 15 -53.48177 -327.61043 513% - 0s
     0 2 -327.57472 0 13 -53.48177 -327.57472 512% - 0s
H 136 322 -199.4028229 -327.57472 64.3% 3.0 0s
H 264 322 -226.5814028 -327.57472 44.6% 3.2 0s
H 1138 806 -300.0734385 -327.57472 9.16% 2.9 0s
H 1578 442 -302.8250805 -326.42894 7.79% 2.9 0s
H33784 1717 -302.8251207 -315.62847 4.23% 2.3 2s
H33791 1636 -311.2578474 -315.62847 1.40% 2.3 2s
 66994 5793 infeasible 91 -311.25785 -315.62847 1.40% 2.4 5s
 184282 9413 infeasible 100 -311.25785 -315.62847 1.40% 2.1 10s
 390579 24443 -315.62847 107 14 -311.25785 -315.62847 1.40% 1.6 15s
 621856 32762 -315.62847 100 9 -311.25785 -315.62847 1.40% 1.4 20s
 820322 32676 infeasible 90 -311.25785 -315.62847 1.40% 1.3 25s
 1054675 32766 -315.62847 95 8 -311.25785 -315.62847 1.40% 1.2 30s
 1271044 32907 infeasible 100 -311.25785 -315.62847 1.40% 1.2 35s
 1483400 33327 -315.62847 115 9 -311.25785 -315.62847 1.40% 1.1 40s
 1691494 32747 -315.62847 96 7 -311.25785 -315.62847 1.40% 1.1 45s
 1876933 33645 infeasible 99 -311.25785 -315.62847 1.40% 1.1 50s
 2030884 33891 infeasible 117 -311.25785 -315.62847 1.40% 1.1 55s

Cutting planes:
  Gomory: 2
  MIR: 38
  StrongCG: 18
  Inf proof: 2

Explored 2206702 nodes (2405517 simplex iterations) in 60.01 seconds (23.12 work units)
Thread count was 8 (of 8 available processors)

Solution count 9: -311.258 -302.825 -302.825 ... -1.36424e-12

Time limit reached
Best objective -3.112578474459e+02, best bound -3.156284748041e+02, gap 1.4042%

User-callback calls 4426781, time in user-callback 1.42 sec

julia> 

julia> Gurobi.GRBwrite(JuMP.unsafe_backend(suspect), "data/suspect.mps")
0

```

After that, I open a new julia REPL and solve it directly, it is very quick.

```julia
julia> import Gurobi

julia> env = Gurobi.Env();
Set parameter Username
Set parameter LicenseID to value 2663413
Academic license - for non-commercial use only - expires 2026-05-11

julia> pModel = Ref{Ptr{Cvoid}}(C_NULL);

julia> Gurobi.GRBreadmodel(env, "data/suspect.mps", pModel)
Read MPS format model from file data/suspect.mps
Reading time = 0.00 seconds
: 28 rows, 108 columns, 207 nonzeros
0

julia> Gurobi.GRBoptimize(pModel[])
Gurobi Optimizer version 12.0.2 build v12.0.2rc0 (win64 - Windows 11.0 (26100.2))

CPU model: 11th Gen Intel(R) Core(TM) i5-1135G7 @ 2.40GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads

Optimize a model with 28 rows, 108 columns and 207 nonzeros
Model fingerprint: 0x8bbf69b7
Variable types: 0 continuous, 108 integer (0 binary)
Coefficient statistics:
  Matrix range [1e+00, 5e+01]
  Objective range [2e+01, 3e+04]
  Bounds range [1e+00, 1e+00]
  RHS range [6e+00, 4e+01]
Presolve removed 7 rows and 7 columns
Presolve time: 0.00s
Presolved: 21 rows, 101 columns, 193 nonzeros
Variable types: 0 continuous, 101 integer (2 binary)
Found heuristic solution: objective 6470.8950980

Root relaxation: objective -3.505153e+02, 53 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 -350.51533 0 11 6470.89510 -350.51533 105% - 0s
H 0 0 1617.7237745 -350.51533 122% - 0s
H 0 0 262.8039698 -350.51533 233% - 0s
H 0 0 46.5482749 -350.51533 853% - 0s
H 0 0 -53.4817681 -350.51533 555% - 0s
     0 0 -339.13976 0 20 -53.48177 -339.13976 534% - 0s
     0 0 -327.61043 0 20 -53.48177 -327.61043 513% - 0s
H 0 0 -99.3727799 -327.61043 230% - 0s
     0 0 -327.61043 0 20 -99.37278 -327.61043 230% - 0s
     0 0 -327.61043 0 19 -99.37278 -327.61043 230% - 0s
H 0 0 -179.8996033 -327.57472 82.1% - 0s
     0 0 -327.57472 0 11 -179.89960 -327.57472 82.1% - 0s
     0 0 -327.57472 0 15 -179.89960 -327.57472 82.1% - 0s
     0 0 -327.57472 0 15 -179.89960 -327.57472 82.1% - 0s
     0 0 -327.57472 0 17 -179.89960 -327.57472 82.1% - 0s
H 0 0 -226.5814028 -327.57472 44.6% - 0s
     0 0 -327.57472 0 25 -226.58140 -327.57472 44.6% - 0s
     0 0 -327.57472 0 17 -226.58140 -327.57472 44.6% - 0s
H 0 0 -251.5349570 -327.57472 30.2% - 0s
     0 0 -327.57472 0 18 -251.53496 -327.57472 30.2% - 0s
     0 0 -327.57472 0 18 -251.53496 -327.57472 30.2% - 0s
     0 0 -327.57472 0 22 -251.53496 -327.57472 30.2% - 0s
     0 0 -327.57472 0 22 -251.53496 -327.57472 30.2% - 0s
     0 2 -327.57472 0 22 -251.53496 -327.57472 30.2% - 0s
H 71 152 -299.8524323 -327.57472 9.25% 1.4 0s
H 1687 920 -303.1600430 -327.57472 8.05% 1.8 0s
H 2084 969 -303.1600841 -324.45429 7.02% 1.8 0s
* 5703 486 104 -307.7487989 -315.62847 2.56% 2.2 0s
*16609 1257 128 -311.2578474 -311.25785 0.00% 2.6 1s

Cutting planes:
  Gomory: 5
  MIR: 15
  StrongCG: 11
  Inf proof: 1

Explored 17226 nodes (44703 simplex iterations) in 1.11 seconds (0.46 work units)
Thread count was 8 (of 8 available processors)

Solution count 10: -311.258 -307.749 -303.16 ... -53.4818
No other solutions better than -311.258

Optimal solution found (tolerance 1.00e-04)
Best objective -3.112578474459e+02, best bound -3.112578474459e+02, gap 0.0000%
0

```

Where is wrong?

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 2, 2025, 12:16pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/2 "2025-07-02T12:16:16Z")

</div>

My first post above itself is an issue—Gurobi’s re-optimize versus optimize-from-scratch.

In this post I reveal an easier issue with Gurobi—it branches a lot and fail to solve a small MIP instance in due time.

The code (independent from my first post) is

```julia
import JuMP, HiGHS
import Gurobi; GRB_ENV = Gurobi.Env();
import LinearAlgebra.⋅ as ⋅
import Random

begin
    seed = 4599192817953253712 # may need to cherry-pick
    Random.seed!(seed)
    J = 9; # number of resources
    W = rand(9:49, J); # the sole width type of each resource
    o = rand(30:0.01:70, J); # occupied cost
    u = 7 * rand(0.75:0.01:1.25, J)/30 .* W; # unit price
    L = rand(6:17, J); # purchasing num limit
    I = 10; # number of customer
    w = rand(7:0.7:24, I); # width needed
    S = 2; # number of scenes
    d = rand(6:35, I, S); # number needed
    Mu = [-65.10070711614117, 15.156225393559112, 0.0, -33369.33855998414, 0.0, 0.0, -32.74174265850863, -114.07026556634082, 256.3700412980869]
    Pi = [617.0164364724739, 306.42333246347755, 188.65466997510032, 462.20679034672725, 539.2412554045152, 539.2412554045152, 208.15789199261513, 316.285736343033, 376.20143138797914]
end

function solve_Rhs(s, Μ, Π, solver)
    if solver == 0
        Rhs = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV));
    else
        Rhs = JuMP.Model(HiGHS.Optimizer);
    end
    JuMP.@variable(Rhs, 0 ≤ b[1:J] ≤ 1, Int);
    JuMP.@variable(Rhs, 0 ≤ m[1:J], Int);
    JuMP.@constraint(Rhs, [j = 1:J], L[j]b[j] ≥ m[j]);
    JuMP.@objective(Rhs, Min, Μ ⋅ b + Π ⋅ m)
    JuMP.@variable(Rhs, 0 ≤ n[1:I, 1:J], Int)
    JuMP.@constraint(Rhs, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
    JuMP.@constraint(Rhs, [j = 1:J], m[j]W[j] ≥ w ⋅ view(n, :, j));
    # JuMP.set_attribute(Rhs, "TimeLimit", 7200)
    JuMP.optimize!(Rhs); JuMP.assert_is_solved_and_feasible(Rhs; allow_local = false)
end

t0 = time()
solve_Rhs(1, Mu, Pi, 1)
println("HiGHSˈtime = $(time() - t0)")

t0 = time()
solve_Rhs(1, Mu, Pi, 0)
println("Gurobiˈtime = $(time() - t0)")

```

On my machine, HiGHS solves this small MIP instance within 11 seconds, while Gurobi cannot solve it within 430 seconds. Here are some results from my machine.

```julia
Solving report
  Status Optimal
  Primal bound -311.257861234
  Dual bound -311.257861234
  Gap 0% (tolerance: 0.01%)
  P-D integral 0.993091429458
  Solution status feasible
                    -311.257861234 (objective)
                    0 (bound viol.)
                    0 (int. viol.)
                    0 (row viol.)
  Timing 9.22 (total)
                    0.00 (presolve)
                    0.00 (solve)
                    0.00 (postsolve)
  Max sub-MIP depth 6
  Nodes 43260
  Repair LPs 0 (0 feasible; 0 iterations)
  LP iterations 100123 (total)
                    742 (strong br.)
                    3653 (separation)
                    14974 (heuristics)

julia> println("HiGHSˈtime = $(time() - t0)")
HiGHSˈtime = 11.315999984741211

julia> solve_Rhs(1, Mu, Pi, 0)
Gurobi Optimizer version 12.0.2 build v12.0.2rc0 (win64 - Windows 11.0 (26100.2))

CPU model: 11th Gen Intel(R) Core(TM) i5-1135G7 @ 2.40GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads

Optimize a model with 28 rows, 108 columns and 207 nonzeros
Model fingerprint: 0xba8bd591
Variable types: 0 continuous, 108 integer (0 binary)
Coefficient statistics:
  Matrix range [1e+00, 5e+01]
  Objective range [2e+01, 3e+04]
  Bounds range [1e+00, 1e+00]
  RHS range [6e+00, 4e+01]
Presolve removed 7 rows and 7 columns
Presolve time: 0.00s
Presolved: 21 rows, 101 columns, 193 nonzeros
Variable types: 0 continuous, 101 integer (2 binary)
Found heuristic solution: objective 6470.8950649

Root relaxation: objective -3.505154e+02, 53 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 -350.51535 0 11 6470.89506 -350.51535 105% - 0s
H 0 0 1617.7237662 -350.51535 122% - 0s
H 0 0 262.8039454 -350.51535 233% - 0s
H 0 0 46.5482567 -350.51535 853% - 0s
H 0 0 -53.4817909 -350.51535 555% - 0s
     0 0 -339.13977 0 20 -53.48179 -339.13977 534% - 0s
     0 0 -327.61045 0 20 -53.48179 -327.61045 513% - 0s
     0 0 -327.61045 0 20 -53.48179 -327.61045 513% - 0s
     0 0 -327.61045 0 19 -53.48179 -327.61045 513% - 0s
H 0 0 -137.1660070 -327.57474 139% - 0s
H 0 0 -161.6096353 -327.57474 103% - 0s
     0 2 -327.57474 0 12 -161.60964 -327.57474 103% - 0s
H 111 172 -199.4028449 -327.57474 64.3% 2.8 0s
H 171 501 -276.4373100 -327.57474 18.5% 2.8 0s
H 1431 734 -303.1600645 -327.57474 8.05% 2.7 0s
H 1440 688 -303.1600726 -327.57474 8.05% 2.7 0s
* 3412 744 118 -307.7488075 -322.12803 4.67% 2.4 0s
* 6870 684 69 -311.2578612 -315.62849 1.40% 2.9 0s
 211849 10731 -315.62849 86 6 -311.25786 -315.62849 1.40% 1.5 5s
 494125 25734 infeasible 88 -311.25786 -315.62849 1.40% 1.4 10s
 759529 36273 -315.62849 97 5 -311.25786 -315.62849 1.40% 1.3 15s
...omit lines.....
 16985734 1127182 -315.62849 88 8 -311.25787 -315.62849 1.40% 1.4 480s
 17121754 1135574 -315.62849 97 11 -311.25787 -315.62849 1.40% 1.4 485s
 17247531 1143422 infeasible 116 -311.25787 -315.62849 1.40% 1.4 490s
...I have no more patience to wait due to the branch size....

```

---

<div class="post-metadata">

**Author:** ![torressa](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/torressa/32/202204_2.png) [@torressa](https://discourse.julialang.org/u/torressa)\
**Post date:** [July 2, 2025, 12:51pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/3 "2025-07-02T12:51:42Z")

</div>

This looks like a case of performance variability [What is performance variability?](https://support.gurobi.com/hc/en-us/articles/26989876720657-What-is-performance-variability)

Solving the same model for three different seeds:

```julia
Explored 64892 nodes (162895 simplex iterations) in 2.05 seconds (0.94 work units)
Explored 89482 nodes (200879 simplex iterations) in 2.94 seconds (1.12 work units)
*I interrupted this* Explored 5586604 nodes (6007635 simplex iterations) in 87.60 seconds (33.64 work units)

```

This can be mitigated by parameter tuning. E.g. setting Presolve=2 gives:

```julia
Explored 15606 nodes (32777 simplex iterations) in 1.04 seconds (0.34 work units)
Explored 7848 nodes (22143 simplex iterations) in 0.46 seconds (0.21 work units)
Explored 32147 nodes (69106 simplex iterations) in 1.44 seconds (0.57 work units)

```

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 2, 2025, 1:41pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/4 "2025-07-02T13:41:58Z")

</div>

Thanks for your opinion. Your info is helpful.

But in my these 2 cases (e.g. the more straightforward second post), the performance is ill acceptable (infinite branches, non-terminating).

> [@torressa](#):
>
> Solving the same model for three different seeds:

Could you provide a JuMP code to set this?

> [@torressa](#):
>
> setting Presolve=2

This indeed works for my second post:

```julia
julia> println("Gurobiˈtime = $(time() - t0)")
Gurobiˈtime = 25.765000104904175

```

But it doesn’t in a general sense help in the algorithm written in my first post. and hence the following question.

Is an MIP like this very hard? Notice there is no Warning in the Logging, i.e. no numeric issues.

```julia
Optimize a model with 28 rows, 108 columns and 207 nonzeros
Model fingerprint: 0xb387f3d4
Variable types: 0 continuous, 108 integer (0 binary)
Coefficient statistics:
  Matrix range [1e+00, 5e+01]
  Objective range [1e+02, 4e+05]
  Bounds range [1e+00, 1e+00]
  RHS range [6e+00, 4e+01]
Presolve removed 6 rows and 6 columns
Presolve time: 0.00s
Presolved: 22 rows, 102 columns, 195 nonzeros
Variable types: 0 continuous, 102 integer (3 binary)

```

The logging for solving it is like this

```julia
  35801363 5038524 -1200.6027 95 8 -1150.0657 -1267.9358 10.2% 2.2 1260s

```

It’s unimaginable for me, until today. Congratulations to me for discovering a class of very hard MIP problems, which I think is instructive. They occurs along the iterations of my algorithm in my first post.

---

<div class="post-metadata">

**Author:** ![torressa](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/torressa/32/202204_2.png) [@torressa](https://discourse.julialang.org/u/torressa)\
**Post date:** [July 2, 2025, 5:08pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/5 "2025-07-02T17:08:21Z")

</div>

You can will have to run the model multiple times changing the Seed parameter and calling reset every time (I can’t find a JuMP function for this but you can always call the `GRBreset` C function from Gurobi.jl)

Presolve=2 was on my machine after a quick experiment, I’m sure better parameters exists. In general, this can happen regardless of the solver, that’s why you should set a termination criteria (MIPGap, TimeLimit).

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 3, 2025, 12:30am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/6 "2025-07-03T00:30:22Z")

</div>

> [@torressa](#):
>
> run the model multiple times changing the Seed parameter and calling reset

I write it as (updated, complete code see post #13)

```julia
rhs_bnd, rhs_val = -Inf, Inf # ✅ Reset memory
for Gurobiˈtry = 1:Gurobiˈtry_TIMES 
    JuMP.optimize!(Rhs)
    JuMP.primal_status(Rhs) == JuMP.FEASIBLE_POINT || error("primal_status ERROR after optimize!")
    termination = JuMP.termination_status(Rhs)
    if termination == JuMP.TIME_LIMIT
        new_rhs_bnd, new_rhs_val = JuMP.objective_bound(Rhs), JuMP.objective_value(Rhs)
        new_rhs_bnd > rhs_bnd && (rhs_bnd = new_rhs_bnd) # update dual side incumbent
        if new_rhs_val < rhs_val # update primal side incumbent
            rhs_val = new_rhs_val
            local_B, local_M = ı.(b), ı.(m) # primal solutions needed
        end
        @info "Fail solving Rhs: $Gurobiˈtry/$Gurobiˈtry_TIMES, re-pick Seed..."
        JuMP.set_attribute(Rhs, "Seed", rand(0:typemax(Int32)))
        err_code = Gurobi.GRBreset(JuMP.unsafe_backend(Rhs), 0)
        err_code == 0 || error("Gurobi reset ERROR: $err_code")
    elseif termination == JuMP.OPTIMAL
        rhs_bnd = JuMP.objective_bound(Rhs)
        local_B, local_M = ı.(b), ı.(m)
        break # [early break]
    else
        error("terminate with $termination")
    end
end

```

it works for my second post case.

> [@torressa](#):
>
> this can happen regardless of the solver

credible. it seems also applies for HiGHS.

(This is quite an advanced usage…)

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 3, 2025, 4:37am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/7 "2025-07-03T04:37:56Z")

</div>

Here is a SmallHardModel having a `6s` TimeLimit, fails to be solved to `OPTIMAL` within 50 consecutive attempts with different Gurobi’s seeds.

What’s your opinion?

> **MPS file of SmallHardModel**
>
> ```julia
> * Signature: 0x34eaa3c19661dfee
> NAME 
> ROWS
> N OBJ
> G R0      
> G R1      
> G R2      
> G R3      
> G R4      
> G R5      
> G R6      
> G R7      
> G R8      
> G R9      
> G R10     
> G R11     
> G R12     
> G R13     
> G R14     
> G R15     
> G R16     
> G R17     
> G R18     
> G R19     
> G R20     
> G R21     
> G R22     
> G R23     
> G R24     
> G R25     
> G R26     
> G R27     
> COLUMNS
> MARKER 'MARKER' 'INTORG'
> b[1] R0 7
> b[2] R1 11
> b[3] OBJ 1.5464237204312122e-02
> b[3] R2 13
> b[4] R3 16
> b[5] R4 11
> b[6] R5 11
> b[7] OBJ -1.3355453911799518e-02
> b[7] R6 9
> b[8] R7 15
> b[9] R8 10
> m[1] OBJ 1.8416499341389969e-01
> m[1] R0 -1
> m[1] R19 45
> m[2] OBJ 0.0934883396545771
> m[2] R1 -1
> m[2] R20 23
> m[3] OBJ 5.6530978622614157e-02
> m[3] R2 -1
> m[3] R21 14
> m[4] OBJ 1.3686853490423276e-01
> m[4] R3 -1
> m[4] R22 34
> m[5] OBJ 1.6307740329014969e-01
> m[5] R4 -1
> m[5] R23 40
> m[6] OBJ 1.6272594593390910e-01
> m[6] R5 -1
> m[6] R24 40
> m[7] OBJ 6.1154026233806148e-02
> m[7] R6 -1
> m[7] R25 15
> m[8] R7 -1
> m[8] R26 23
> m[9] OBJ 1.1317008683069976e-01
> m[9] R8 -1
> m[9] R27 28
> n[1,1] R9 1
> n[1,1] R19 -8.4
> n[2,1] R10 1
> n[2,1] R19 -12.6
> n[3,1] R11 1
> n[3,1] R19 -12.6
> n[4,1] R12 1
> n[4,1] R19 -9.1
> n[5,1] R13 1
> n[5,1] R19 -17.5
> n[6,1] R14 1
> n[6,1] R19 -16.1
> n[7,1] R15 1
> n[7,1] R19 -15.4
> n[8,1] R16 1
> n[8,1] R19 -21.7
> n[9,1] R17 1
> n[9,1] R19 -13.3
> n[10,1] R18 1
> n[10,1] R19 -9.1
> n[1,2] R9 1
> n[1,2] R20 -8.4
> n[2,2] R10 1
> n[2,2] R20 -12.6
> n[3,2] R11 1
> n[3,2] R20 -12.6
> n[4,2] R12 1
> n[4,2] R20 -9.1
> n[5,2] R13 1
> n[5,2] R20 -17.5
> n[6,2] R14 1
> n[6,2] R20 -16.1
> n[7,2] R15 1
> n[7,2] R20 -15.4
> n[8,2] R16 1
> n[8,2] R20 -21.7
> n[9,2] R17 1
> n[9,2] R20 -13.3
> n[10,2] R18 1
> n[10,2] R20 -9.1
> n[1,3] R9 1
> n[1,3] R21 -8.4
> n[2,3] R10 1
> n[2,3] R21 -12.6
> n[3,3] R11 1
> n[3,3] R21 -12.6
> n[4,3] R12 1
> n[4,3] R21 -9.1
> n[5,3] R13 1
> n[5,3] R21 -17.5
> n[6,3] R14 1
> n[6,3] R21 -16.1
> n[7,3] R15 1
> n[7,3] R21 -15.4
> n[8,3] R16 1
> n[8,3] R21 -21.7
> n[9,3] R17 1
> n[9,3] R21 -13.3
> n[10,3] R18 1
> n[10,3] R21 -9.1
> n[1,4] R9 1
> n[1,4] R22 -8.4
> n[2,4] R10 1
> n[2,4] R22 -12.6
> n[3,4] R11 1
> n[3,4] R22 -12.6
> n[4,4] R12 1
> n[4,4] R22 -9.1
> n[5,4] R13 1
> n[5,4] R22 -17.5
> n[6,4] R14 1
> n[6,4] R22 -16.1
> n[7,4] R15 1
> n[7,4] R22 -15.4
> n[8,4] R16 1
> n[8,4] R22 -21.7
> n[9,4] R17 1
> n[9,4] R22 -13.3
> n[10,4] R18 1
> n[10,4] R22 -9.1
> n[1,5] R9 1
> n[1,5] R23 -8.4
> n[2,5] R10 1
> n[2,5] R23 -12.6
> n[3,5] R11 1
> n[3,5] R23 -12.6
> n[4,5] R12 1
> n[4,5] R23 -9.1
> n[5,5] R13 1
> n[5,5] R23 -17.5
> n[6,5] R14 1
> n[6,5] R23 -16.1
> n[7,5] R15 1
> n[7,5] R23 -15.4
> n[8,5] R16 1
> n[8,5] R23 -21.7
> n[9,5] R17 1
> n[9,5] R23 -13.3
> n[10,5] R18 1
> n[10,5] R23 -9.1
> n[1,6] R9 1
> n[1,6] R24 -8.4
> n[2,6] R10 1
> n[2,6] R24 -12.6
> n[3,6] R11 1
> n[3,6] R24 -12.6
> n[4,6] R12 1
> n[4,6] R24 -9.1
> n[5,6] R13 1
> n[5,6] R24 -17.5
> n[6,6] R14 1
> n[6,6] R24 -16.1
> n[7,6] R15 1
> n[7,6] R24 -15.4
> n[8,6] R16 1
> n[8,6] R24 -21.7
> n[9,6] R17 1
> n[9,6] R24 -13.3
> n[10,6] R18 1
> n[10,6] R24 -9.1
> n[1,7] R9 1
> n[1,7] R25 -8.4
> n[2,7] R10 1
> n[2,7] R25 -12.6
> n[3,7] R11 1
> n[3,7] R25 -12.6
> n[4,7] R12 1
> n[4,7] R25 -9.1
> n[5,7] R13 1
> n[5,7] R25 -17.5
> n[6,7] R14 1
> n[6,7] R25 -16.1
> n[7,7] R15 1
> n[7,7] R25 -15.4
> n[8,7] R16 1
> n[8,7] R25 -21.7
> n[9,7] R17 1
> n[9,7] R25 -13.3
> n[10,7] R18 1
> n[10,7] R25 -9.1
> n[1,8] R9 1
> n[1,8] R26 -8.4
> n[2,8] R10 1
> n[2,8] R26 -12.6
> n[3,8] R11 1
> n[3,8] R26 -12.6
> n[4,8] R12 1
> n[4,8] R26 -9.1
> n[5,8] R13 1
> n[5,8] R26 -17.5
> n[6,8] R14 1
> n[6,8] R26 -16.1
> n[7,8] R15 1
> n[7,8] R26 -15.4
> n[8,8] R16 1
> n[8,8] R26 -21.7
> n[9,8] R17 1
> n[9,8] R26 -13.3
> n[10,8] R18 1
> n[10,8] R26 -9.1
> n[1,9] R9 1
> n[1,9] R27 -8.4
> n[2,9] R10 1
> n[2,9] R27 -12.6
> n[3,9] R11 1
> n[3,9] R27 -12.6
> n[4,9] R12 1
> n[4,9] R27 -9.1
> n[5,9] R13 1
> n[5,9] R27 -17.5
> n[6,9] R14 1
> n[6,9] R27 -16.1
> n[7,9] R15 1
> n[7,9] R27 -15.4
> n[8,9] R16 1
> n[8,9] R27 -21.7
> n[9,9] R17 1
> n[9,9] R27 -13.3
> n[10,9] R18 1
> n[10,9] R27 -9.1
> MARKER 'MARKER' 'INTEND'
> RHS
> RHS1 R9 20
> RHS1 R10 16
> RHS1 R11 21
> RHS1 R12 22
> RHS1 R13 6
> RHS1 R14 8
> RHS1 R15 21
> RHS1 R16 19
> RHS1 R17 35
> RHS1 R18 20
> BOUNDS
> UP BND1 b[1] 1
> UP BND1 b[2] 1
> UP BND1 b[3] 1
> UP BND1 b[4] 1
> UP BND1 b[5] 1
> UP BND1 b[6] 1
> UP BND1 b[7] 1
> UP BND1 b[8] 1
> UP BND1 b[9] 1
> LI BND1 m[1] 0
> LI BND1 m[2] 0
> LI BND1 m[3] 0
> LI BND1 m[4] 0
> LI BND1 m[5] 0
> LI BND1 m[6] 0
> LI BND1 m[7] 0
> LI BND1 m[8] 0
> LI BND1 m[9] 0
> LI BND1 n[1,1] 0
> LI BND1 n[2,1] 0
> LI BND1 n[3,1] 0
> LI BND1 n[4,1] 0
> LI BND1 n[5,1] 0
> LI BND1 n[6,1] 0
> LI BND1 n[7,1] 0
> LI BND1 n[8,1] 0
> LI BND1 n[9,1] 0
> LI BND1 n[10,1] 0
> LI BND1 n[1,2] 0
> LI BND1 n[2,2] 0
> LI BND1 n[3,2] 0
> LI BND1 n[4,2] 0
> LI BND1 n[5,2] 0
> LI BND1 n[6,2] 0
> LI BND1 n[7,2] 0
> LI BND1 n[8,2] 0
> LI BND1 n[9,2] 0
> LI BND1 n[10,2] 0
> LI BND1 n[1,3] 0
> LI BND1 n[2,3] 0
> LI BND1 n[3,3] 0
> LI BND1 n[4,3] 0
> LI BND1 n[5,3] 0
> LI BND1 n[6,3] 0
> LI BND1 n[7,3] 0
> LI BND1 n[8,3] 0
> LI BND1 n[9,3] 0
> LI BND1 n[10,3] 0
> LI BND1 n[1,4] 0
> LI BND1 n[2,4] 0
> LI BND1 n[3,4] 0
> LI BND1 n[4,4] 0
> LI BND1 n[5,4] 0
> LI BND1 n[6,4] 0
> LI BND1 n[7,4] 0
> LI BND1 n[8,4] 0
> LI BND1 n[9,4] 0
> LI BND1 n[10,4] 0
> LI BND1 n[1,5] 0
> LI BND1 n[2,5] 0
> LI BND1 n[3,5] 0
> LI BND1 n[4,5] 0
> LI BND1 n[5,5] 0
> LI BND1 n[6,5] 0
> LI BND1 n[7,5] 0
> LI BND1 n[8,5] 0
> LI BND1 n[9,5] 0
> LI BND1 n[10,5] 0
> LI BND1 n[1,6] 0
> LI BND1 n[2,6] 0
> LI BND1 n[3,6] 0
> LI BND1 n[4,6] 0
> LI BND1 n[5,6] 0
> LI BND1 n[6,6] 0
> LI BND1 n[7,6] 0
> LI BND1 n[8,6] 0
> LI BND1 n[9,6] 0
> LI BND1 n[10,6] 0
> LI BND1 n[1,7] 0
> LI BND1 n[2,7] 0
> LI BND1 n[3,7] 0
> LI BND1 n[4,7] 0
> LI BND1 n[5,7] 0
> LI BND1 n[6,7] 0
> LI BND1 n[7,7] 0
> LI BND1 n[8,7] 0
> LI BND1 n[9,7] 0
> LI BND1 n[10,7] 0
> LI BND1 n[1,8] 0
> LI BND1 n[2,8] 0
> LI BND1 n[3,8] 0
> LI BND1 n[4,8] 0
> LI BND1 n[5,8] 0
> LI BND1 n[6,8] 0
> LI BND1 n[7,8] 0
> LI BND1 n[8,8] 0
> LI BND1 n[9,8] 0
> LI BND1 n[10,8] 0
> LI BND1 n[1,9] 0
> LI BND1 n[2,9] 0
> LI BND1 n[3,9] 0
> LI BND1 n[4,9] 0
> LI BND1 n[5,9] 0
> LI BND1 n[6,9] 0
> LI BND1 n[7,9] 0
> LI BND1 n[8,9] 0
> LI BND1 n[9,9] 0
> LI BND1 n[10,9] 0
> ENDATA
> 
> ```

On my machine it looks like

```julia
julia> Gurobi.GRBreadmodel(env, "data/SmallHardModel.mps", pModel)
Read MPS format model from file data/SmallHardModel.mps
Reading time = 0.00 seconds
: 28 rows, 108 columns, 207 nonzeros
0

julia> Gurobi.GRBoptimize(pModel[])
Gurobi Optimizer version 12.0.2 build v12.0.2rc0 (win64 - Windows 11.0 (26100.2))

CPU model: 11th Gen Intel(R) Core(TM) i5-1135G7 @ 2.40GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads

Optimize a model with 28 rows, 108 columns and 207 nonzeros
Model fingerprint: 0x7f6ed3f7
Variable types: 0 continuous, 108 integer (0 binary)
Coefficient statistics:
  Matrix range [1e+00, 5e+01]
  Objective range [1e-02, 2e-01]
  Bounds range [1e+00, 1e+00]
  RHS range [6e+00, 4e+01]
Presolve removed 8 rows and 9 columns
Presolve time: 0.00s
Presolved: 20 rows, 99 columns, 190 nonzeros
Variable types: 0 continuous, 99 integer (1 binary)
Found heuristic solution: objective 8.7316200

Root relaxation: objective 8.532225e+00, 39 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 8.53223 0 11 8.73162 8.53223 2.28% - 0s
H 0 0 8.7124249 8.53223 2.07% - 0s
H 0 0 8.6558940 8.53223 1.43% - 0s
H 0 0 8.6550018 8.53268 1.41% - 0s
     0 0 8.53584 0 21 8.65500 8.53584 1.38% - 0s
H 0 0 8.5938478 8.53590 0.67% - 0s
H 31 67 8.5624172 8.53635 0.30% 3.5 0s
H 1774 1205 8.5537813 8.53668 0.20% 2.4 0s
H 1775 1146 8.5509696 8.53779 0.15% 2.4 0s
 13068305 147141 8.54043 136 8 8.54218 8.54043 0.02% 1.6 325s
omit lines...

```

In this solve, it can only converges to 0.02% after 2s, but the gap never reduces further. Will you recommend me to broaden the MIPGap to above 0.02%?

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 3, 2025, 4:49am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/8 "2025-07-03T04:49:41Z")

</div>

> [@WalterMadelim](#):
>
> Will you recommend me to broaden the MIPGap to above 0.02%?

Well, thought it doesn’t look decent. I do succeed finishing my Benders’ algorithm algorithm 😅.

> **Code of my algorithm**
>
> ```julia
> import JuMP.value as ı
> import LinearAlgebra.⋅ as ⋅
> import Random, JuMP, Gurobi; GRB_ENV = Gurobi.Env();
> 
> seed = 4599192817953253712 # may need to cherry-pick
> Random.seed!(seed)
> J = 9; # number of resources
> W = rand(9:49, J); # the sole width type of each resource
> o = rand(30:0.01:70, J); # occupied cost
> u = 7 * rand(0.75:0.01:1.25, J)/30 .* W; # unit price
> L = rand(6:17, J); # purchasing num limit
> I = 10; # number of customer
> w = rand(7:0.7:24, I); # width needed
> S = 2; # number of scenes
> d = rand(6:35, I, S); # number needed
> 
> if false # This is the original MIP. I attempt to use Benders Decomposition to solve it
> Xtn = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # extensive formulation
> # Xtn = JuMP.Model(HiGHS.Optimizer)
> JuMP.@variable(Xtn, 0 ≤ b[1:J] ≤ 1, Int); # use this site or not
> JuMP.@variable(Xtn, 0 ≤ m[1:J], Int); # order an amount
> JuMP.@constraint(Xtn, [j = 1:J], L[j]b[j] ≥ m[j]); # source num is limited
> JuMP.@expression(Xtn, Obj_1, o ⋅ b + u ⋅ m);
> # each s indexes a block
> JuMP.@variable(Xtn, 0 ≤ m2[1:J, 1:S]);
> JuMP.@variable(Xtn, 0 ≤ n[1:I, 1:J, 1:S], Int);
> JuMP.@constraint(Xtn, [i = 1:I, s = 1:S], sum(view(n, i, :, s)) ≥ d[i, s]);
> JuMP.@constraint(Xtn, [j = 1:J, s = 1:S], m2[j, s]W[j] ≥ w ⋅ view(n, :, j, s)); 
> JuMP.@constraint(Xtn, [j = 1:J, s = 1:S], m2[j, s] == m[j]); # 🟡 Linking
> # temporarily don't consider obj from 2nd stage
> JuMP.@objective(Xtn, Min, Obj_1 + sum(0/S for s = 1:S)) 
> JuMP.unset_integer.(b)
> JuMP.unset_integer.(m)
> JuMP.unset_integer.(n)
> JuMP.optimize!(Xtn); JuMP.assert_is_solved_and_feasible(Xtn; allow_local = false)
> LP_val = JuMP.objective_value(Xtn)
> JuMP.set_integer.(b)
> JuMP.set_integer.(m)
> JuMP.set_integer.(n)
> JuMP.optimize!(Xtn); JuMP.assert_is_solved_and_feasible(Xtn; allow_local = false)
> IP_val = JuMP.objective_value(Xtn);
> @info "$LP_val < $IP_val" # 949.2709333333335 < 986.3943333333333
> @info "optimal solution" B = JuMP.value.(b) M = JuMP.value.(m)
> # ┌ Info: optimal solution
> # │ B =
> # │ 9-element Vector{Float64}:
> # │ 1.0
> # │ 1.0
> # │ 0.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ M =
> # │ 9-element Vector{Float64}:
> # │ 7.0
> # │ 11.0
> # │ 0.0
> # │ 16.0
> # │ 8.0
> # │ 11.0
> # │ 9.0
> # │ 15.0
> # └ 10.0
> end;
> 
> function Min_s_prim_model(s) # evaluate the 2nd-stage-scene-s exactly ⚠️ Feasible Region is NOT fixed
> Min_s = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV))
> JuMP.@variable(Min_s, m2[1:J]) # JuMP.fix.(m2, m) 🟡 [Jit] Linking 
> JuMP.@variable(Min_s, 0 ≤ n[1:I, 1:J], Int)
> JuMP.@constraint(Min_s, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
> JuMP.@constraint(Min_s, [j = 1:J], m2[j]W[j] ≥ w ⋅ view(n, :, j))
> JuMP.@objective(Min_s, Min, 0/S)
> return Min_s
> end; function Rhs_model(s) # 🟢 Feasible Region is fixed
> Rhs = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV))
> JuMP.set_attribute(Rhs, "TimeLimit", 6)
> JuMP.set_attribute(Rhs, "MIPGap", 1e-3)
> JuMP.@variable(Rhs, 0 ≤ b[1:J] ≤ 1, Int)
> JuMP.@variable(Rhs, 0 ≤ m[1:J], Int); JuMP.@constraint(Rhs, [j = 1:J], L[j]b[j] ≥ m[j])
> JuMP.@objective(Rhs, Min, 0) # TODO JuMP.set_objective_coefficient
> # For the single scene s
> JuMP.@variable(Rhs, 0 ≤ n[1:I, 1:J], Int)
> JuMP.@constraint(Rhs, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
> JuMP.@constraint(Rhs, [j = 1:J], m[j]W[j] ≥ w ⋅ view(n, :, j))
> return Rhs
> end; function Rhs_surrogate_plus_s(s)
> Sur = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV));
> JuMP.@variable(Sur, μ[1:J])
> JuMP.@variable(Sur, π[1:J])
> JuMP.@constraint(Sur, [true; μ; π] ∈ JuMP.MOI.NormOneCone(1+2J))
> JuMP.@variable(Sur, θ) # > Rhs(μ, π) initially
> JuMP.@objective(Sur, Max, θ) # need more 🟣 θ-Cuts (can be reused)
> JuMP.set_objective_coefficient.(Sur, μ, -falses(J)) # -B
> JuMP.set_objective_coefficient.(Sur, π, -falses(J)) # -M
> return Sur
> end; macro ready_Rhs_and_Sur_code()
> esc(quote
> Sur_model_vec = [Rhs_surrogate_plus_s(s) for s = 1:S]; # ✅ used to generate trial coefficients of Lag cut
> Rhs_model_vec = [Rhs_model(s) for s = 1:S]; # ✅ bring the RHS of Lag's cut
> JuMP.set_silent.(Sur_model_vec)
> JuMP.set_silent.(Rhs_model_vec)
> for s = 1:S # Add an initial round of cut to ensure ⚫ BOUNDEDNESS
> Μ = falses(J); Μ[1] = true; Π = falses(J) # 💡 some valid trial multipliers
> local b, m, B, M, Rhs
> Rhs = Rhs_model_vec[s]
> b, m = Rhs[:b], Rhs[:m]
> JuMP.set_objective_coefficient.(Rhs, b, Μ); JuMP.set_objective_coefficient.(Rhs, m, Π)
> JuMP.optimize!(Rhs); JuMP.assert_is_solved_and_feasible(Rhs; allow_local = false)
> B, M = JuMP.value.(b), JuMP.value.(m)
> Sur = Sur_model_vec[s]
> # because this is the very 1st cut, we obviate the need for COT-test
> JuMP.@constraint(Sur, Sur[:θ] ≤ B ⋅ Sur[:μ] + M ⋅ Sur[:π]) # 🟣 θ-Cuts
> JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false) # ⚫ BOUNDEDNESS
> end
> end)
> end
> 
> @ready_Rhs_and_Sur_code()
> prim_model_vec = [Min_s_prim_model(s) for s = 1:S]; # ✅ evaluate the 2nd stage ObjVals
> JuMP.set_silent.(prim_model_vec);
> 
> Mst = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # ✅ used to generate trial 1st-stage solutions (only Int solutions admit valid cut-offs)
> JuMP.@variable(Mst, 0 ≤ b[1:J] ≤ 1, Int);
> JuMP.@variable(Mst, 0 ≤ m[1:J], Int); JuMP.@constraint(Mst, [j = 1:J], L[j]b[j] ≥ m[j]);
> JuMP.@objective(Mst, Min, o ⋅ b + u ⋅ m);
> JuMP.set_silent(Mst);
> JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
> B, M = ı.(b), ı.(m); # to be globally ✂️ off
> SmallHardModel = nothing
> while true # iteratively find a violating cut for the Master problem
> infeasible_block = 0
> for s = 1:S # each block
> Min_s = prim_model_vec[s]; JuMP.fix.(Min_s[:m2], M)
> JuMP.optimize!(Min_s); termination = JuMP.termination_status(Min_s)
> if termination == JuMP.INFEASIBLE
> infeasible_block = s
> break # [early break]
> elseif termination != JuMP.OPTIMAL
> error("primal 2nd stage: $termination")
> end
> end
> s = infeasible_block
> if s == 0
> @info "😃 An OPTIMAL solution is derived!" ObjVal=(o ⋅ B + u ⋅ M) B M
> break
> end
> s == 0 && error("We have now derived ")
> Sur = Sur_model_vec[s]; θ, μ, π = Sur[:θ], Sur[:μ], Sur[:π];
> JuMP.set_objective_coefficient.(Sur, μ, -B)
> JuMP.set_objective_coefficient.(Sur, π, -M)
> JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false); ub = JuMP.objective_bound(Sur)
> Μ, Π = ı.(μ), ı.(π); # to be locally ✂️ off
> Θ = ı(θ); # [🪜Upper] to be locally ✂️ off
> rhs = -Inf
> while true # try generate the coefficient of Lag cut until success
> local b, m # This is NOT the (b, m) of the outermost Master
> Rhs = Rhs_model_vec[s]; b, m = Rhs[:b], Rhs[:m]
> JuMP.set_objective_coefficient.(Rhs, b, Μ); JuMP.set_objective_coefficient.(Rhs, m, Π)
> optimal = false
> for ø = 1:50
> JuMP.optimize!(Rhs)
> termination = JuMP.termination_status(Rhs)
> if termination == JuMP.TIME_LIMIT
> @info "Fail solving Rhs: $ø/100, re-pick Seed..."
> JuMP.set_attribute(Rhs, "Seed", rand(0:typemax(Int32)))
> err_code = Gurobi.GRBreset(JuMP.unsafe_backend(Rhs), 0)
> err_code == 0 || error("Gurobi reset ERROR: $err_code")
> elseif termination == JuMP.OPTIMAL
> optimal = true
> break # [early break]
> else
> error("terminate with $termination")
> end
> end
> if optimal == false
> @error "The Rhs as SmallHardModel is too hard to be solved"
> SmallHardModel = Rhs
> return SmallHardModel
> end
> rhs = JuMP.objective_value(Rhs)
> lb = rhs - B ⋅ Μ - M ⋅ Π # here (B, M) are from the outermost Master
> if lb > 4e-3 # 🟣 global COT
> @info "Find a cut for master with vio = $lb"
> break # go back to add cut for the outer Master
> end
> local_B, local_M = ı.(b), ı.(m)
> vertical_distance = Θ - (local_B ⋅ Μ + local_M ⋅ Π) # 🪜 having this since you approximate `Rhs_model_s`
> if vertical_distance > 8e-4 # 🟣 local COT
> JuMP.@constraint(Sur, θ ≤ local_B ⋅ μ + local_M ⋅ π)
> else
> error("weird 568435686586442")
> end
> JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false); ub = JuMP.objective_bound(Sur)
> Μ, Π = ı.(μ), ı.(π); # to be locally ✂️ off
> Θ = ı(θ); # [🪜Upper] to be locally ✂️ off
> end
> JuMP.@constraint(Mst, Μ ⋅ b + Π ⋅ m ≥ rhs) # 🟣 global COT
> JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
> B, M = ı.(b), ı.(m); # to be globally ✂️ off
> gub = o ⋅ B + u ⋅ M
> @info "Master ub = $gub"
> end
> 
> ```

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [July 3, 2025, 5:00am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/9 "2025-07-03T05:00:43Z")

</div>

I think you have a fundamental misconception about solving MIPs.

They are NP hard, but we have algorithms that are surprisingly effective in practice.

Nonetheless, there are no performance guarantees. A model may solve in 1 second, but, with a small modification, never solve in 1000 years. Solving with various seeds is one trick to exploit this performance variability, but it does not mean that there exists a seed which will solve fast.

I do **not** recommend your time-limit loop approach with the seeds. David’s solution was a way to demonstrate the performance variability of MIPs, it was not an algorithm. You might equally be better off using the same number of seconds with a single solve.

Re the MIP gap: this is problem and application dependent. It is up to you to decide whether it matter for your problem. We can’t provide general advice.

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 3, 2025, 5:22am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/10 "2025-07-03T05:22:40Z")

</div>

> [@odow](#):
>
> have algorithms that are surprisingly effective in practice.

I think the majority is cutting planes, secondarily it’s a set of primal heuristics.  
Branching is the last fallback—they finally become hopeless.

> [@odow](#):
>
> there are no performance guarantees

I was just wondering “Ah, you Gurobi cannot solve a model with 20 rows, 99 columns, 190 nonzeros”??? It was unimaginable for me.

> [@odow](#):
>
> You might equally be better off using the same number of seconds with a single solve.

I didn’t get you here. What is “using the same number of seconds”? I’m already doing this.

> [@odow](#):
>
> I do **not** recommend your time-limit loop approach with the seeds.

I find it works somehow. I update my successful code in my #7 post. You can check it out. On my computer it runs about 5 minutes and the global optimal solution is found. (This is a Benders decomposition for 2 stage stochastic programming having general integer variables on both stages).

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [July 3, 2025, 5:31am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/11 "2025-07-03T05:31:13Z")

</div>

> [@WalterMadelim](#):
>
> Ah, you Gurobi cannot solve a model with 20 rows, 99 columns, 190 nonzeros”??? It was unimaginable for me.

Yes. This should be expected. There are NO performance guarantees, even for tiny MIPs.

---

<div class="post-metadata">

**Author:** ![torressa](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/torressa/32/202204_2.png) [@torressa](https://discourse.julialang.org/u/torressa)\
**Post date:** [July 3, 2025, 1:05pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/12 "2025-07-03T13:05:30Z")

</div>

For example, many open instances in the MIPLIB have a modest number of variables/constraints: [The Open Tag](https://miplib.zib.de/tag_open.html). Similar for the [hard instances](https://miplib.zib.de/tag_hard.html).

It sometimes helps to try to reformulate the original problem with an alternative formulation in the hopes that it is easier to solve

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 3, 2025, 3:55pm UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/13 "2025-07-03T15:55:24Z")

</div>

Thanks for sharing the links.

> [@torressa](#):
>
> reformulate the original problem with an alternative formulation in the hopes that it is easier to solve

that entails very special domain knowledge, I think.

> [@odow](#):
>
> I do **not** recommend your time-limit loop approach with the seeds. David’s solution was a way to demonstrate the performance variability of MIPs, it was not an algorithm. You might equally be better off using the same number of seconds with a single solve.

I now have a deeper understanding about this.

Firstly, it is indeed helpful to re-try, **even if** you try **only twice**. You may likely to fail if you just hinges on a single solve. But the situation might be greatly improved if `Seed + GRBreset` is used. I’ll update my time-limit loop code—I wrote a better one. You can revisit my code therein. As an initial trial, we can set `Gurobiˈtry_TIMES = 1`, which is equivalent to an ordinary single solve. If the results are not desired, we may have good luck by promoting `Gurobiˈtry_TIMES` to `2`, `3`…

My previous “successful” code involved setting an attribute `MIPGap` to `1e-3` (10 times of the Gurobi’s default). Since I don’t think it is some praiseworthy trick, I managed to invent a new version in which the **only 2** non-default attributes are `"Seed"` and `"TimeLimit"`.

Although it takes longer to solve to global optimality, I believe this version of implementation will be numerically way more stable.

> **Code without deteriorating MIPGap**
>
> ```julia
> import JuMP.value as ı
> import LinearAlgebra.⋅ as ⋅
> import Random, JuMP, Gurobi; GRB_ENV = Gurobi.Env();
> 
> seed = 4599192817953253712 # may need to cherry-pick
> Random.seed!(seed)
> J = 9; # number of resources
> W = rand(9:49, J); # the sole width type of each resource
> o = rand(30:0.01:70, J); # occupied cost
> u = 7 * rand(0.75:0.01:1.25, J)/30 .* W; # unit price
> L = rand(6:17, J); # purchasing num limit
> I = 10; # number of customer
> w = rand(7:0.7:24, I); # width needed
> S = 2; # number of scenes
> d = rand(6:35, I, S); # number needed
> 
> if false # This is the original MIP. I attempt to use Benders Decomposition to solve it
> Xtn = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # extensive formulation
> JuMP.@variable(Xtn, 0 ≤ b[1:J] ≤ 1, Int); # use this site or not
> JuMP.@variable(Xtn, 0 ≤ m[1:J], Int); # order an amount
> JuMP.@constraint(Xtn, [j = 1:J], L[j]b[j] ≥ m[j]);
> JuMP.@expression(Xtn, Obj_1, o ⋅ b + u ⋅ m);
> # each s indexes a block
> JuMP.@variable(Xtn, 0 ≤ m2[1:J, 1:S]);
> JuMP.@variable(Xtn, 0 ≤ n[1:I, 1:J, 1:S], Int); # 2nd-stage dispatch variable
> JuMP.@constraint(Xtn, [i = 1:I, s = 1:S], sum(view(n, i, :, s)) ≥ d[i, s]);
> JuMP.@constraint(Xtn, [j = 1:J, s = 1:S], m2[j, s]W[j] ≥ w ⋅ view(n, :, j, s));
> JuMP.@constraint(Xtn, [j = 1:J, s = 1:S], m2[j, s] == m[j]); # 🟡 Linking
> JuMP.@objective(Xtn, Min, Obj_1 + sum(0/S for s = 1:S)) 
> JuMP.unset_integer.(b)
> JuMP.unset_integer.(m)
> JuMP.unset_integer.(n)
> JuMP.optimize!(Xtn); JuMP.assert_is_solved_and_feasible(Xtn; allow_local = false)
> LP_val = JuMP.objective_value(Xtn)
> JuMP.set_integer.(b)
> JuMP.set_integer.(m)
> JuMP.set_integer.(n)
> JuMP.optimize!(Xtn); JuMP.assert_is_solved_and_feasible(Xtn; allow_local = false)
> IP_val = JuMP.objective_value(Xtn);
> @info "LP_ObjVal = $LP_val < $IP_val = IP_ObjVal" # 949.2709333333335 < 986.3943333333333
> @info "optimal solution" B = JuMP.value.(b) M = JuMP.value.(m)
> # ┌ Info: optimal solution
> # │ B =
> # │ 9-element Vector{Float64}:
> # │ 1.0
> # │ 1.0
> # │ 0.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ 1.0
> # │ M =
> # │ 9-element Vector{Float64}:
> # │ 7.0
> # │ 11.0
> # │ 0.0
> # │ 16.0
> # │ 8.0
> # │ 11.0
> # │ 9.0
> # │ 15.0
> # └ 10.0
> end; function Min_s_prim_model(s) # evaluate the 2nd-stage-scene-s exactly ⚠️ Feasible Region is NOT fixed
> Min_s = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV))
> JuMP.@variable(Min_s, m2[1:J]) # JuMP.fix.(m2, m) 🟡 [Jit] Linking 
> JuMP.@variable(Min_s, 0 ≤ n[1:I, 1:J], Int)
> JuMP.@constraint(Min_s, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
> JuMP.@constraint(Min_s, [j = 1:J], m2[j]W[j] ≥ w ⋅ view(n, :, j))
> JuMP.@objective(Min_s, Min, 0/S)
> return Min_s
> end; function Rhs_model(s) # 🟢 Feasible Region is fixed
> # [■ Pivot] since:
> # 1. To update Sur, we need feasible solution, feasible ObjVal
> # 2. To update Master, we need ObjBnd
> Rhs = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV))
> JuMP.set_attribute(Rhs, "TimeLimit", 6)
> JuMP.@variable(Rhs, 0 ≤ b[1:J] ≤ 1, Int)
> JuMP.@variable(Rhs, 0 ≤ m[1:J], Int); JuMP.@constraint(Rhs, [j = 1:J], L[j]b[j] ≥ m[j])
> JuMP.@objective(Rhs, Min, 0) # TODO JuMP.set_objective_coefficient
> # For the single scene s
> JuMP.@variable(Rhs, 0 ≤ n[1:I, 1:J], Int)
> JuMP.@constraint(Rhs, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s])
> JuMP.@constraint(Rhs, [j = 1:J], m[j]W[j] ≥ w ⋅ view(n, :, j))
> return Rhs
> end; function Rhs_surrogate_plus_s(s)
> Sur = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV));
> JuMP.@variable(Sur, μ[1:J])
> JuMP.@variable(Sur, π[1:J])
> JuMP.@constraint(Sur, [true; μ; π] ∈ JuMP.MOI.NormOneCone(1+2J))
> JuMP.@variable(Sur, θ) # > Rhs(μ, π) initially
> JuMP.@objective(Sur, Max, θ) # need more 🟣 θ-Cuts (can be reused)
> JuMP.set_objective_coefficient.(Sur, μ, -falses(J)) # -B
> JuMP.set_objective_coefficient.(Sur, π, -falses(J)) # -M
> return Sur
> end; macro check_feas_return_s()
> esc(quote
> infeasible_block = 0 
> for s = 1:S # each block
> Min_s = prim_model_vec[s]; JuMP.fix.(Min_s[:m2], M)
> JuMP.optimize!(Min_s); termination = JuMP.termination_status(Min_s)
> if termination == JuMP.INFEASIBLE
> infeasible_block = s
> break # [early break]
> elseif termination != JuMP.OPTIMAL
> error("primal 2nd stage: $termination")
> end
> end; 
> s = infeasible_block
> if s == 0
> @info "😃 A global optimal solution is derived!" ObjVal=(o ⋅ B + u ⋅ M) B M
> return Mst # the well-trained Master
> end
> end)
> end; macro optimize_Sur_return_Mu_Pi_Theta()
> esc(quote
> JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false)
> Μ, Π = ı.(μ), ı.(π); # to be locally ✂️ off
> Θ = ı(θ); # [🪜Upper] to be locally ✂️ off
> end)
> end
> 
> Sur_model_vec = [Rhs_surrogate_plus_s(s) for s = 1:S]; # ✅ used to generate trial coefficients of Lag cut
> Rhs_model_vec = [Rhs_model(s) for s = 1:S]; # ✅ bring the RHS of Lag's cut
> JuMP.set_silent.(Sur_model_vec)
> JuMP.set_silent.(Rhs_model_vec)
> for s = 1:S # Add an initial round of cut to ensure ⚫ BOUNDEDNESS
> Μ = falses(J); Μ[1] = true; Π = falses(J) # 💡 some valid trial multipliers
> local b, m, Rhs
> Rhs = Rhs_model_vec[s]
> b, m = Rhs[:b], Rhs[:m]
> JuMP.set_objective_coefficient.(Rhs, b, Μ); JuMP.set_objective_coefficient.(Rhs, m, Π)
> JuMP.optimize!(Rhs); JuMP.assert_is_solved_and_feasible(Rhs; allow_local = false)
> local_B, local_M = JuMP.value.(b), JuMP.value.(m)
> Sur = Sur_model_vec[s]
> # because this is the very 1st cut, we obviate the need for COT-test
> JuMP.@constraint(Sur, Sur[:θ] ≤ local_B ⋅ Sur[:μ] + local_M ⋅ Sur[:π]) # 🟣 θ-Cuts
> JuMP.optimize!(Sur); JuMP.assert_is_solved_and_feasible(Sur; allow_local = false) # ⚫ BOUNDEDNESS
> end
> prim_model_vec = [Min_s_prim_model(s) for s = 1:S]; # ✅ evaluate the 2nd stage ObjVals
> JuMP.set_silent.(prim_model_vec);
> 
> Random.seed!(658354866) # 🔯 change this seed brings different stochastic process during solving
> Mst = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # ✅ used to generate trial 1st-stage solutions (only Int solutions admit valid cut-offs)
> JuMP.@variable(Mst, 0 ≤ b[1:J] ≤ 1, Int);
> JuMP.@variable(Mst, 0 ≤ m[1:J], Int); JuMP.@constraint(Mst, [j = 1:J], L[j]b[j] ≥ m[j]);
> JuMP.@objective(Mst, Min, o ⋅ b + u ⋅ m);
> JuMP.set_silent(Mst);
> JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
> B, M = ı.(b), ı.(m); # to be globally ✂️ off
> rhs_bnd, local_B, local_M = -Inf, zeros(J), zeros(J) # memory
> Gurobiˈtry_TIMES = 1 # 🔯
> global_COT = 5.0
> while true # iteratively find a violating cut for the Master problem
> @check_feas_return_s()
> Sur = Sur_model_vec[s]; θ, μ, π = Sur[:θ], Sur[:μ], Sur[:π];
> Rhs = Rhs_model_vec[s]; local_b, local_m = Rhs[:b], Rhs[:m]
> JuMP.set_objective_coefficient.(Sur, μ, -B)
> JuMP.set_objective_coefficient.(Sur, π, -M)
> @optimize_Sur_return_Mu_Pi_Theta()
> local_COT, local_saturation = 2.0, false # # ✅ Reset memory
> while true # try generate the coefficient of Lag cut until success
> JuMP.set_objective_coefficient.(Rhs, local_b, Μ)
> JuMP.set_objective_coefficient.(Rhs, local_m, Π)
> rhs_bnd, rhs_val = -Inf, Inf # ✅ Reset memory
> for Gurobiˈtry = 1:Gurobiˈtry_TIMES # 🔯
> JuMP.optimize!(Rhs)
> JuMP.primal_status(Rhs) == JuMP.FEASIBLE_POINT || error("primal_status 5437473754")
> termination = JuMP.termination_status(Rhs)
> if termination == JuMP.TIME_LIMIT
> new_rhs_bnd, new_rhs_val = JuMP.objective_bound(Rhs), JuMP.objective_value(Rhs)
> new_rhs_bnd > rhs_bnd && (rhs_bnd = new_rhs_bnd) # update dual side incumbent
> if new_rhs_val < rhs_val # update primal side incumbent
> rhs_val = new_rhs_val
> local_B, local_M = ı.(local_b), ı.(local_m) # primal solutions needed
> end
> @info "Fail solving Rhs: $Gurobiˈtry/$Gurobiˈtry_TIMES, re-pick Seed..."
> JuMP.set_attribute(Rhs, "Seed", rand(0:typemax(Int32)))
> err_code = Gurobi.GRBreset(JuMP.unsafe_backend(Rhs), 0)
> err_code == 0 || error("Gurobi reset ERROR: $err_code")
> elseif termination == JuMP.OPTIMAL
> rhs_bnd = JuMP.objective_bound(Rhs)
> local_B, local_M = ı.(local_b), ı.(local_m)
> break # [early break]
> else
> error("terminate with $termination")
> end
> end
> lb = rhs_bnd - B ⋅ Μ - M ⋅ Π # here (B, M) are from the outermost Master
> if lb > global_COT # 🟣 
> @info "Find a cut for master with vio = $lb"
> break
> elseif local_saturation
> if global_COT ≥ 1e-4
> global_COT /= 2
> else
> error("There is no more violating cuts, in ALL models")
> end
> end
> ### ### ### ### ### ### ### ### ### ### ### ### ### 
> vertical_distance = Θ - (local_B ⋅ Μ + local_M ⋅ Π) # 🪜 having this since you approximate `Rhs_model_s`
> local_saturation = false
> if vertical_distance > local_COT # 🟣 
> JuMP.@constraint(Sur, θ ≤ local_B ⋅ μ + local_M ⋅ π)
> elseif local_COT ≥ 1e-4
> local_COT /= 2
> else
> local_saturation = true
> end
> @optimize_Sur_return_Mu_Pi_Theta()
> end
> JuMP.@constraint(Mst, Μ ⋅ b + Π ⋅ m ≥ rhs_bnd)
> JuMP.optimize!(Mst); JuMP.assert_is_solved_and_feasible(Mst; allow_local = false)
> B, M = ı.(b), ı.(m); # to be globally ✂️ off
> @info "Master ub = $(o ⋅ B + u ⋅ M)"
> end
> 
> ```

---

<div class="post-metadata">

**Author:** ![WalterMadelim](https://avatars.discourse-cdn.com/v4/letter/w/3e96dc/32.png) [@WalterMadelim](https://discourse.julialang.org/u/WalterMadelim)\
**Post date:** [July 6, 2025, 5:34am UTC](https://discourse.julialang.org/t/gurobi-gets-stuck-after-loaded-mip-start-from-previous-solve/130406/14 "2025-07-06T05:34:33Z")

</div>

> [@WalterMadelim](#):
>
> Here is a SmallHardModel having a `6s` TimeLimit, fails to be solved to `OPTIMAL` within 50 consecutive attempts with different Gurobi’s seeds.
> 
> What’s your opinion?

From my experience gathered these days, I seriously suspect that the difficulty of the MIP instance is due in part to the data.  
The data like `1.5464237204312122e-02` appears to have **too many** digits.  
Plus the integer restrictions, the default termination rgap `0.01%` might be overly strict.

But data like this are not a physical or economical coefficient from the **original** formulation, but rather from the dual variable (i.e., generated by the cut generation algorithm along the iterations). Therefore, writing an algorithm ourselves is bound to be more intricate than directly employing a black-box solver.

* * *

Here I offer a toy example, in the following code (`🍅` line), if `digits = 3`, exact optimality can be reached within `2` seconds. If you modify it to `digits = 4`, then Gurobi get stuck at `0.02%` rgap and branches endlessly. (To add one additional information, it seems HiGHS perform better sometimes)

```julia
import LinearAlgebra.⋅ as ⋅
import JuMP, Gurobi; GRB_ENV = Gurobi.Env();
import Random

seed = 5026110663637712176
Random.seed!(seed)
J = 9; # number of resources
W = rand(9:49, J); # the sole width type of each resource
o = rand(30:0.01:70, J); # occupied cost
u = 7 * rand(0.75:0.01:1.25, J)/30 .* W; # unit price
L = rand(6:17, J); # purchasing num limit
I = 10; # number of customer
w = rand(7:0.7:24, I); # width needed
S = 2; # number of scenes
d = rand(6:35, I, S); # number needed
c = rand(0:0.001:0.1, I, J); # transport cost

(s, Πb, Πm, Π0, an_UB) = (2, [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.15629861471878964, 0.32235412781500744, 0.06320282559299903, 0.0, 0.10588629037089686, 0.0, 0.12587639701895473, 0.22638174448335244], 0.0, 6.448866244173365)
Πm = round.(Πm; digits = 3) # 🍅 digits = 4 is much harder

Rhs = JuMP.Model(() -> Gurobi.Optimizer(GRB_ENV)); # 7137655 344699 cutoff 123 6.44887 6.44772 0.02% 2.5 145s
JuMP.set_attribute(Rhs, "MIPGap", 0);
JuMP.set_attribute(Rhs, "MIPGapAbs", 0);
JuMP.@variable(Rhs, 0 ≤ b[1:J] ≤ 1, Int);
JuMP.@variable(Rhs, 0 ≤ m[1:J], Int); JuMP.@constraint(Rhs, [j = 1:J], L[j]b[j] ≥ m[j]);
JuMP.@variable(Rhs, 0 ≤ n[1:I, 1:J], Int);
JuMP.@objective(Rhs, Min, +(Πb ⋅ b, Πm ⋅ m, Π0/S * c ⋅ n));
JuMP.@constraint(Rhs, [j = 1:J], m[j]W[j] ≥ w ⋅ view(n, :, j));
JuMP.@constraint(Rhs, [i = 1:I], sum(view(n, i, :)) ≥ d[i, s]);
JuMP.optimize!(Rhs);
JuMP.MOI.get(JuMP.unsafe_backend(Rhs), Gurobi.ModelAttribute("MIPGap"))

```

@odow I don’;t know what is `JuMP.unsafe_backend`, but I find it usable.