Ah. If all the b_dict
are 1
, then you have a few options, and you could think up a variety of reformulations.
With Gurobi, you could use the GeometricMeanCone: Standard form · MathOptInterface
@variable(model, t)
x = [p[r] * ratio[r] for r in od]
@constraint(model, vcat(t, x) in MOI.GeometricMeanCone(length(od) + 1))
@objective(model, Max, t)
but when I ran it, I encountered numerical issues, likely because of the big-M constant.
You could also take logs and solve the MINLP directly, using Juniper:
using JuMP, Juniper, Ipopt
model = Model(
optimizer_with_attributes(
Juniper.Optimizer,
"nl_solver" => optimizer_with_attributes(Ipopt.Optimizer),
),
)
# ...
@NLobjective(model, Max, sum(log(p[r] * ratio[r]) for r in od))
You could also try other things, like a lazy constraint cutting plane method:
@variable(model, t[od] <= M)
@objective(model, Max, sum(t))
function my_callback(cb_data)
for r in od
t_r = callback_value(cb_data, t[r])
ratio_r = callback_value(cb_data, ratio[r])
x = p[r] * ratio_r
f = log(x)
if t_r > f + 1e-6
con = @build_constraint(t[r] <= f + 1 / x * (ratio[r] - ratio_r))
MOI.submit(model, MOI.LazyConstraint(cb_data), con)
return
end
end
return
end
MOI.set(model, MOI.LazyConstraintCallback(), my_callback)
optimize!(model)