JuMP solution is not correct: PiecewiseLinearOpt

Specific Domains Optimization (Mathematical)
1

Good morning,

I’m working on a MILP model where a non-linear function is linearized using PiecewiseLinearOpt and Gurobi.
The model is solved and the objective value reached is correct, but the constraints are not met.
For example:
constraint Demand[1]: g_gen[1] + g_gen_h[1] = 351
but the solution is g_gen_h[1] = 0 ; g_gen[1] = 101

JuMP.value(g_gen_h[1]) is not getting the correct value.

The code is available on https://github.com/salvaguerrero/NonConvexPLF_code/blob/master/hydro.jl

model = Model()
set_optimizer(model, Gurobi.Optimizer)
@variable(model, g_gen[t in Time]   >= 0)
@variable(model, g_gen_h[t in Time] >= 0)
@variable(model, w_tur[t in Time]   >= 0)
@variable(model, w_lev[t in 0:T_num]>= 0)
@variable(model, w_spi[t in Time]   >= 0)

@objective(model, Min, sum(g_gen[t]*Opex for t in Time) )

@constraint(model,Demand[t in Time], g_gen_h[t] + g_gen[t] == dem[t] ) 

@constraint(model, [t in Time], L_lim <= g_gen[t] <= U_lim)
@constraint(model, [t in Time], L_lim <= g_gen_h[t] <= U_lim)

@constraint(model, Reservoir[t in Time],  w_lev[t] == w_lev[t-1] + w_apo -(w_tur[t] + w_spi[t]) )

@constraint(model, w_lev[0] == w_lev_0 )
@constraint(model, w_lev[T_num] == w_lev_f )

p(q,v) = 9.81*q*v#/a/b*eff

q_range = w_tur_m:1:w_tur_M #0:128
v_range = w_lev_m:1:w_lev_M #0:128
for t in Time
   g_gen_h[t]  = piecewiselinear(model, w_tur[t] , w_lev[t], v_range, q_range, p)
end

optimize!(model)

Thank you in advance, Salvador

2

piecewiselinear returns a variable. In the for-loop, you are over-writing it, rather than adding a constraint that they must be equal.

Do instead something like this:

using JuMP
import Gurobi
import PiecewiseLinearOpt

Opex  = 25
L_lim, U_lim = 0, 250
Time = 1:24
dem = [200 * sin(t * 2 * pi / 24) + 300 for t in Time] 
w_lev_0, w_lev_f = 10, 10
w_apo = 20
w_lev_m, w_lev_M, w_tur_m, w_tur_M = 0, 128, 0, 128
p(q, v) = 9.81 * q * v
q_range = w_tur_m:w_tur_M
v_range = w_lev_m:w_lev_M

model = Model(Gurobi.Optimizer)
@variables(model, begin
    L_lim  <= g_gen[Time] <= U_lim
    w_tur[Time] >= 0
    w_lev[vcat(0, Time)] >= 0
    w_spi[Time] >= 0
end)
fix(w_lev[0], w_lev_0; force = true)
fix(w_lev[end], w_lev_f; force = true)
g_gen_h = [
    PiecewiseLinearOpt.piecewiselinear(model, w_tur[t] , w_lev[t], v_range, q_range, p)
    for t in Time
]
set_lower_bound.(g_gen_h, L_lim)
set_upper_bound.(g_gen_h, U_lim)
@objective(model, Min, Opex * sum(g_gen[t] for t in Time))
@constraints(model, begin
    [t in Time], g_gen_h[t] + g_gen[t] == dem[t]
    [t in Time],  w_lev[t] == w_lev[t-1] + w_apo - (w_tur[t] + w_spi[t])
end) 
optimize!(model)
1 Like