Efficient approach to modifying existing JuMP constraints

Hello everyone,

Long time lurker, first time poster. :smiley: I’m trying to add a new expression to an existing JuMP constraint, but my approach seems pretty slow. Am I approaching this wrong? Thank you in advance.

# @expression(m, New_Expr[i in set[:key_1]], Ξ± * Another_Expr[i])

for i in set[:key_1], t in set[:time_24h]
    # ...store the original contraint function, (=> slow execution)
    constraint = constraint_object(m[:Existing_Cstr][i, t])
    # ...delete it from the model,
    delete(m, m[:Existing_Cstr][i, t])
    # ...and reinsert the modified constraint
    m[:Existing_Cstr][i, t] = @constraint(m, constraint.func - New_Expr[i] == 0.)
end
2 Likes

What are the constraint and expression types generally speaking? Affine/linear, quadratic and/or nonlinear?

1 Like

Hi @aartinian, welcome :smile:

This depends one what, specifically, you want to modify. Some examples:

In general, it is possible to write something that is much more efficient than deleting and re-adding a constraint.

To give some specific advice, can you provide a reproducible example? It doesn’t have to be your exact code; something similar will do.

1 Like

Thank you both for your lightning-fast response to my question. I’ll prepare a reproducible example and post it later today.

1 Like

Parameters are another good way to modify constraints GitHub - jump-dev/ParametricOptInterface.jl: Extension for dealing with parameters

2 Likes

While preparing the example below, I realised that I was actually using direct_model to initialise the model. This was necessary due to some additional adjustments of Gurobi. However, it is no longer needed. Switching to the Model constructor removes the bottleneck mentioned above. Still, I would be very interested in a more efficient solution here.

using JuMP, Gurobi

function main()

    sto = TimerOutput()

    var_1 = 1
    var_2 = 1:100
    var_3 = 1:10
    var_4 = 1:144
    var_5 = [:a]
    var_6 = [1, 2, 3]

    sets = Dict(
        :A => [var_1],
        :B => Dict(var_1 => collect(var_2)),
        :C => Dict((i => collect(var_3) for i in var_2)...),
        :D => collect(var_4),
        :E => var_5,
        :F => Dict(var_1 => [:b1, :b2, :b3]),
        :G => Dict((j => var_6 for j in var_3)...),
        :H => Dict((j => 100 + 10 * j) for j in var_3),
        :I => Dict(((j, k) => 0.1 * k + 0.05 * j) for j in var_3 for k in var_6)
    )

    # model = Model(() ->  Gurobi.Optimizer(GRB_ENV))
    model = direct_model(Gurobi.Optimizer(GRB_ENV))

    # variables
    @variable(model, 0 <= X[sets[:A], i in sets[:B][var_1], j in sets[:C][i], t in sets[:D]] <= sets[:H][j])
    @variable(model, 0 <= Y[sets[:A], n in sets[:E], br in sets[:F][var_1]] <= 100)

    # expressions
    @expression(model, Z[a in sets[:A], i in sets[:B][var_1], j in sets[:C][i], b in sets[:G][j], t in sets[:D]], X[a, i, j, t] * sets[:I][(j, b)])
    @expression(model, W[a in sets[:A], n in sets[:E]], sum(Y[a, n, br] for br in sets[:F][var_1]))

    # constraints
    @constraint(model, U[a in sets[:A], b in var_6, t in sets[:D]], sum(Z[a, i, j, b, t] for i in sets[:B][var_1], j in sets[:C][i] if b in sets[:G][j]) == 0)

    # ... a few moments later ...

    calculated_parameter = 0.5
    @expression(model, V[a in sets[:A], n in sets[:E]], calculated_parameter * W[a, n])

    for a in sets[:A], t in sets[:D]
        @timeit sto "get" constraint = constraint_object(model[:U][a, 1, t])
        @timeit sto "del" delete(model, model[:U][a, 1, t])
        @timeit sto "set" model[:U][a, 1, t] = @constraint(
            model, constraint.func - V[a, :a] == 0.
        )
    end

    print_timer(sto)
end

main()

With direct_model :

 ────────────────────────────────────────────────────────────────────
                            Time                    Allocations      
                   ───────────────────────   ────────────────────────
 Tot / % measured:      42.1s /  93.1%           1.40GiB /   2.3%    

 Section   ncalls     time    %tot     avg     alloc    %tot      avg
 ────────────────────────────────────────────────────────────────────
 get          144    39.1s   99.7%   272ms   16.2MiB   49.9%   115KiB
 set          144    105ms    0.3%   731ΞΌs   16.2MiB   50.0%   115KiB
 del          144   16.1ms    0.0%   112ΞΌs   6.75KiB    0.0%    48.0B
 ────────────────────────────────────────────────────────────────────

With the Model constructor:

  ────────────────────────────────────────────────────────────────────
                            Time                    Allocations      
                   ───────────────────────   ────────────────────────
 Tot / % measured:      2.97s /   0.8%           1.40GiB /   2.0%    

 Section   ncalls     time    %tot     avg     alloc    %tot      avg
 ────────────────────────────────────────────────────────────────────
 get          144   11.9ms   52.2%  82.4ΞΌs   12.3MiB   42.2%  87.2KiB
 set          144   10.7ms   47.0%  74.3ΞΌs   16.7MiB   57.6%   119KiB
 del          144    182ΞΌs    0.8%  1.26ΞΌs   38.7KiB    0.1%     275B
 ────────────────────────────────────────────────────────────────────

What version of Gurobi.jl do you have installed? Try updating to the latest version, we recently fixed a performance issue that is probably causing this.

I would probably do something like this:

    # @expression(model, W[a in sets[:A], n in sets[:E]], sum(Y[a, n, br] for br in sets[:F][var_1]))
    @variable(model, W[a in sets[:A], n in sets[:E]])
    @constraint(
        model, 
        [a in sets[:A], n in sets[:E]],
        W[a, n] == sum(Y[a, n, br] for br in sets[:F][var_1]),
    )
    @constraint(
        model, 
        U[a in sets[:A], b in var_6, t in sets[:D]],
        sum(Z[a, i, j, b, t] for i in sets[:B][var_1], j in sets[:C][i] if b in sets[:G][j]) == 0,
    )
    calculated_parameter = 0.5
    for a in sets[:A], t in sets[:D]
        set_normalized_coefficient(model[:U][a, 1, t], W[a, :a], -calculated_parameter)
    end

where you introduce a new variable to represent the expression, and then you modify the coefficient of that variable in the relevant constraints.

Thank you, @odow . It works perfectly and even opens the door for further improvements. :smile:

Currently using Gurobi.jl 1.0.0, but will definitely update soon as part of switching from Julia 1.8 to 1.11.

You don’t have to switch versions of Julia to update a package. Just run:

julia> import Pkg

julia> Pkg.update("Gurobi")
    Updating registry at `~/.julia/registries/General.toml`
    Updating `~/.julia/dev/Convex/grb/Project.toml`
  [2e9cd046] ↑ Gurobi v1.0.0 β‡’ v1.2.3
    Updating `~/.julia/dev/Convex/grb/Manifest.toml`
  [2e9cd046] ↑ Gurobi v1.0.0 β‡’ v1.2.3

There are no breaking changes, just performance and bug fixes:

2 Likes