(Jump, Gurobi)Binary Variables end up with values no equal 1 or 0

Hirnfolter · July 31, 2021, 6:01pm

Hello,
i try to opimise something and have the variable

    @variable(m, Taskstartstation[1:Stationnumber, 1:Ordernumber, 1:Tasknumbertotal], binary=true);

or

@variable(m, Taskstartstation[1:Stationnumber, 1:Ordernumber, 1:Tasknumbertotal], bin);

and somehow the solution is for example.

julia>     Oldtaskstartstation = JuMP.value.(Taskstartstation)
6×10×13 Array{Float64, 3}:
[:, :, 1] =
  0.0  0.0         0.0  0.0         0.0          0.0  0.0         0.0         0.0         0.0
  0.0  0.0         0.0  0.0         0.0          0.0  0.0         0.0         0.0         0.0
  1.0  3.37834e-6  0.0  7.99834e-6  7.99834e-6   1.0  0.999997    7.99834e-6  7.99834e-6  7.99834e-6
  0.0  0.0         0.0  0.0         0.0          0.0  0.0         0.0         0.0         0.0
  0.0  0.0         0.0  0.0         0.0          0.0  0.0         0.0         0.0         0.0
 -0.0  0.999997    1.0  0.999992    0.999992    -0.0  3.37834e-6  0.999992    0.999992    0.999992

Why are there values not equal to 1 or 0?!
Thanks in advance.

DrChainsaw · July 31, 2021, 9:00pm

Someone else can probably give a more accurate answer, but I think the underlying solver considers the result within tolerances. If Jump reports the problem as solved you should be able to safely round the result. As you can see, the numbers are very close to either 0 or 1.

Hirnfolter · July 31, 2021, 9:11pm

@DrChainsaw

Thanks for the fast answer.
The thing is, that it considered 0,99999 equals 1 was confusing.
If somone finds a solution how to stop it, it would be great.

I ran in problems because if this, later in the progamm.
Because i used simple high values like 99999*bin var in some constraints and now when the variable is not 0 but 0.000001, the constraint gives a false true.

So a solution would be good. My workaround was taking the 99999 and reduced it to 99, so the tolerance fits, but its shady.

DrChainsaw · July 31, 2021, 9:24pm

Are the constraints violated even when you round the values to integers?

If not, I’m 78.9999992% certain that the result is good (I’m not very experienced with mathematical programming hence the low confidence).

This issue might also give further insight: https://github.com/jump-dev/Cbc.jl/issues/141

Edit: Link from issue which seems extra relevant: https://www.gurobi.com/documentation/9.0/refman/num_grb_tolerances_and_the.html

odow · July 31, 2021, 10:01pm

Hi there! Since this is your first post, welcome.

Why are there values not equal to 1 or 0?!

Solvers find solutions to a certain tolerance. They don’t compute using exact arithmetic. Your solution looks a little bad, typically you would expect 1e-8 instead of 1e-6, but we’d need the full model to explain why. Do you have the output of the Gurobi log?

Since it seems you might be new to numerical optimization, Gurobi has a good documentation section on this and many related topics. I highly recommend the reading:

Hirnfolter · August 4, 2021, 5:53pm

Can i simply take a integer and limit it to =>0 and <=1

cvanaret · August 4, 2021, 6:26pm

Could be that Gurobi didn’t branch (it creates subproblems in which a binary variable is set to 0 and 1, respectively) but applied numerical techniques (such as cutting planes) that are subject to numerical errors.

odow · August 4, 2021, 7:05pm

No. That will have the same issue.

Can you provide the log or a minimal working example?

Hirnfolter · August 5, 2021, 3:01pm

thats the problem mainly

Slackhighnumber = 9999

    Slacklownumber = 0.000001

    Stationnumber = size(Stationtask,1)

    Tasknumbertotal = size(Taskmodel[1],2)

    Tasknumber = Tasknumbertotal - 2

    Ordernumber = size(Ordersubset,1)

    Start = Tasknumbertotal

    End = Tasknumbertotal - 1

    MaxAGV = AVGnumber

    Transporttime = 0

    AGVEndtoStart = 0 #Time it takes the AGV to get from the last Station to the first Station

    Oldordersnumber = size(Oldtaskstart,1)



    m = nothing

    m = Model(optimizer_with_attributes(Solver))

    set_optimizer_attribute(m, "TimeLimit", Timelimitins)

    #Variable

    @variable(m, Endtime >= 0, Int);

    @variable(m, Taskstart[1:Ordernumber, 1:Tasknumbertotal] >= 0, Int); #Divided everywhere by /100 as a Precisionlimit

    @variable(m, Taskstartstation[1:Stationnumber, 1:Ordernumber, 1:Tasknumbertotal], Bin);

    @variable(m, T1beforeT2[1:Ordernumber, 1:Tasknumbertotal, 1:Tasknumbertotal], Bin);

    @variable(m, OT1beforeOT2[1:Ordernumber,1:Ordernumber, 1:Tasknumbertotal, 1:Tasknumbertotal], Bin);

    @variable(m, OrderstartedbeforeOrder[1:Ordernumber,1:Ordernumber], Bin);

    @variable(m, OrdernotfinshedbeforerOrder[1:Ordernumber,1:Ordernumber], Bin);

    @variable(m, ActiveOrder[1:Ordernumber,1:Ordernumber], Bin);

    @variable(m, Transport[1:Ordernumber,1:Tasknumbertotal,1:Tasknumbertotal], Bin);



    #Fix

    

    for Index in CartesianIndices(Oldtaskstart)

        fix(Taskstart[Index[1], Index[2]], Oldtaskstart[Index[1], Index[2]]*100; force = true)

    end 

    for Index in CartesianIndices(Oldtaskstartstation)

        fix(Taskstartstation[Index[1], Index[2], Index[3]], Oldtaskstartstation[Index[1], Index[2], Index[3]]; force = true)

    end 

    

    #complete all Orders ASAP

    @objective(m, Min, Endtime + sum(Taskstart[O, End] for O = 1:Ordernumber)*Slacklownumber);

    #Objective sub Constraint to find Maximum Endtime

    @constraint(m, Endtime1[O=1:Ordernumber,T=1:Tasknumber], Endtime >= 0.01 * Taskstart[O, End]);



    #Taskoverlapp + Transporttime

    @constraint(m, Taskoverlap1[O=1:Ordernumber,T1=1:(Tasknumbertotal-1),T2=(T1+1):Tasknumbertotal], 0.01 * Taskstart[O,T1]+Taskmodel[Ordersubset[O]][T1]+Transport[O,T1,T2]*Transporttime<=0.01 * Taskstart[O,T2]+T1beforeT2[O,T1,T2]*Slackhighnumber);

    @constraint(m, Taskoverlap2[O=1:Ordernumber,T1=1:(Tasknumbertotal-1),T2=(T1+1):Tasknumbertotal], 0.01 * Taskstart[O,T2]+Taskmodel[Ordersubset[O]][T2]+Transport[O,T1,T2]*Transporttime<=0.01 * Taskstart[O,T1]+(T1beforeT2[O,T1,T2]-1)*-Slackhighnumber);

    @constraint(m, Transport1[S=1:Stationnumber,O=1:Ordernumber,T1=1:(Tasknumbertotal-1),T2=(T1+1):Tasknumbertotal],Taskstartstation[S,O,T1]-Taskstartstation[S,O,T2]<=Transport[O,T1,T2]);

    

    #Constraint Start and End

    @constraint(m, TaskstartStart[O=1:Ordernumber,T=1:Tasknumber], 0.01 * Taskstart[O, Start] <= 0.01 * Taskstart[O, T]);

    @constraint(m, TaskstartEnd[O=1:Ordernumber,T=1:Tasknumber], 0.01 * Taskstart[O, T] + Taskmodel[Ordersubset[O]][T] <= 0.01 * Taskstart[O, End]);



    #Constraint Precedence

    @constraint(m, Precedencecon[O=1:Ordernumber,P=1:size(Precedence[Ordersubset[O]],1)], 0.01 * Taskstart[O, Precedence[Ordersubset[O]][P][1]] <= 0.01 * Taskstart[O, Precedence[Ordersubset[O]][P][2]]);



    #Constraint every Starttask has 1 Station

    @constraint(m, TaskStation1[O=1:Ordernumber,T=1:Tasknumbertotal],sum(Taskstartstation[S,O,T] for S=1:Stationnumber) == 1);



    #Constraint Stationoverlap  ////Start is not on the checked Station = (Taskstartstation[S,O1,T1]-1)*(-1)*Slackhighnumber+(Taskstartstation[S,O2,T2]-1)*(-1)*Slackhighnumber) 

    @constraint(m, Stationoverlap1[S=1:Stationnumber,O1=1:(Ordernumber-1), O2=(O1+1):Ordernumber,T1=1:Tasknumbertotal,T2=1:Tasknumbertotal], 

    0.01 * Taskstart[O1,T1]+Taskmodel[Ordersubset[O1]][T1]<=0.01 * Taskstart[O2,T2]+OT1beforeOT2[O1,O2,T1,T2]*Slackhighnumber + (Taskstartstation[S,O1,T1]-1)*-Slackhighnumber+(Taskstartstation[S,O2,T2]-1)*-Slackhighnumber);

    @constraint(m, Stationoverlap2[S=1:Stationnumber,O1=1:(Ordernumber-1), O2=(O1+1):Ordernumber,T1=1:Tasknumbertotal,T2=1:Tasknumbertotal], 

    0.01 * Taskstart[O2,T2]+Taskmodel[Ordersubset[O2]][T2]<=0.01 * Taskstart[O1,T1]+(OT1beforeOT2[O1,O2,T1,T2]-1)*-Slackhighnumber + (Taskstartstation[S,O1,T1]-1)*-Slackhighnumber+(Taskstartstation[S,O2,T2]-1)*-Slackhighnumber);

    

    #Constraint Starttask at corrosponding Station

    @constraint(m, StationwithTask[S=1:Stationnumber, O=1:Ordernumber,T=1:Tasknumbertotal],Taskstartstation[S,O,T] <= Stationtask[S][T]);



    #Constraint maximum avgs

    @constraint(m, OrderstartedbeforeOrder1[O1=1:Ordernumber,O2 in 1:Ordernumber; O2!=O1],0.01 * Taskstart[O1,Start]-0.01 * Taskstart[O2,Start] + Slacklownumber <= OrderstartedbeforeOrder[O1,O2]*Slackhighnumber);

    @constraint(m, OrdernotfinshedbeforerOrder1[O1=1:Ordernumber,O2=1:Ordernumber; O2!=O1],0.01 * Taskstart[O2,End]-0.01 * Taskstart[O1,Start] + AGVEndtoStart <= OrdernotfinshedbeforerOrder[O1,O2]*Slackhighnumber);

    @constraint(m, ActiveOrder1[O1=1:Ordernumber,O2=1:Ordernumber; O2!=O1],ActiveOrder[O1,O2] >= OrdernotfinshedbeforerOrder[O1,O2] + OrderstartedbeforeOrder[O1,O2] - 1);

    @constraint(m, MaxAVGs1[O1=1:Ordernumber],sum(ActiveOrder[O1,O2] for O2 in 1:Ordernumber if O2!=O1)<=MaxAGV-1);



    #constraint path

    @constraint(m, MovinginX[S1=1:(Stationnumber-1),S2=(S1+1):(Stationnumber),T=1:Tasknumber],Taskstartstation[S1,:,T] + Taskstartstation[S2,:,T+1].<= Stationpos[S2][1] - Stationpos[S1][1] +2);



    #Constraint Taskstart after old Orders

    @constraint(m, OldOrders[O=(Oldordersnumber+1):Ordernumber], 0.01 * Taskstart[O, Start] >= FinisholdOrders);



    optimize!(m);

    

    println("Termination status: ", JuMP.termination_status(m));

    println("Primal status: ", JuMP.primal_status(m));

    println("Objective value: ", JuMP.objective_value(m));

the binary variable are to determine if a process is for or after another process