JuMP - Update constraint with matrix coefficient

Specific Domains Optimization (Mathematical)
2

Hi @luccafaro, welcome to the forum :smile:

A few comments.

About set_normalized_coefficient

The first is that set_normalized_coefficient(con, alpha, H) doesn’t work because con is a Vector of constraints. In general, if you encounter a MethodError in Julia because you are trying to apply a scalar operation over a vector:

julia> set_normalized_coefficient(con, alpha, H)
ERROR: MethodError: no method matching set_normalized_coefficient(::Vector{ConstraintRef{…}}, ::Vector{VariableRef}, ::Matrix{Float64})

you need to use Julia’s broadcasting:

set_normalized_coefficient.(con, alpha, H)

Note the .

However, this isn’t quite correct, because of the normalized part of the function name. JuMP will normalize all constraints by moving the variable terms to the left-hand side. Thus, your normalized constraint is really:

@constraint(model, con, [Up; Uf] .- H * alpha .== 0)

so you need to do

set_normalized_coefficient.(con, alpha, -H)

However, this call is a bit complicated, because you’re broadcasting over a vector of constraints (con) and a vector of variables (alpha) and a matrix of numbers (H).

Because you need con down the rows of H and alpha across the columns of H, you need to transpose the alpha vector so that it is a row. That gives:

set_normalized_coefficient.(con, alpha', -H)

Using broadcasting with a mix of shapes can be tricky to get correct. Note that Julia is not MATLAB, so loops are not slower than vectorized code. You could just write instead:

for i in eachindex(con)
    for j in eachindex(alpha)
        set_normalized_coefficient(con[i], alpha[j], -H[i, j])
    end
end

Using parameters

However, there’s a much easier (and in this case, more efficient way): use Nonlinear Modeling · JuMP

model = Model(Ipopt.Optimizer)
@variable(model, Uf[1:32])
@variable(model, alpha[1:40])
@variable(model, H_param[i in 1:40, j in 1:40] in Parameter(H[i, j]))
@constraint(model, con, [Up; Uf] .== H_param * alpha)
optimize!(model)
@assert is_solved_and_feasible(model)
set_parameter_value.(H_param, 2 * H)
optimize!(model)
@assert is_solved_and_feasible(model)
4 Likes