# Max of a vector as an objective for JuMP?

**URL:** https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911
**Category:** Optimization (Mathematical)
**Tags:** jump
**Created:** [July 11, 2024, 1:49pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911 "2024-07-11T13:49:54Z")
**Posts on this page:** 12
**Page:** 1

<div class="post-metadata">

### Author: ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)
#### Post date: [July 11, 2024, 1:49pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/1 "2024-07-11T13:49:54Z")

</div>

Hi Everyone,

I have an optimization problem similar to [this one](https://github.com/epirecipes/sir-julia/blob/master/markdown/function_map_lockdown_jump/function_map_lockdown_jump.md), except my objective is now the maximum element in a vector (which I want to minimize). Can JuMP handle this as an objective, and if so, what’s the best way to implement this?

---

<div class="post-metadata">

### Author: ![barucden](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/barucden/32/26154_2.png) [@barucden](https://discourse.julialang.org/u/barucden)
#### Post date: [July 11, 2024, 2:48pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/2 "2024-07-11T14:48:40Z")

</div>

Could you introduce an auxiliary variable, say c, add a set of conditions x\_i \leq c for each variable x\_i, and set the objective as \min c?

---

<div class="post-metadata">

### Author: ![blegat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/blegat/32/217090_2.png) [@blegat](https://discourse.julialang.org/u/blegat)
#### Post date: [July 11, 2024, 3:47pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/3 "2024-07-11T15:47:38Z")

</div>

We have special support for it with the [`MOI.NormInfinityCone`](https://jump.dev/JuMP.jl/stable/moi/reference/standard_form/#MathOptInterface.NormInfinityCone) so you can do

```julia
@variable(model, c)
@constraint(model, [c; x] in MOI.NormInfinityCone(length(x) + 1))
@objective(model, Min, c)

```

This will be reformulated to a formulation equivalent to what @barucden suggests except if you use the solver [Hypatia](https://github.com/jump-dev/Hypatia.jl/) which supports this cone natively.

---

<div class="post-metadata">

### Author: ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)
#### Post date: [July 11, 2024, 5:55pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/4 "2024-07-11T17:55:21Z")

</div>

Thanks to both @barucden and @blegat! I also have a constraint where I have a variable `@variable(model, υ[1:(T+1)])`, where I want to have a constraint where the count of variables \> 0 is less than some maximum value. I tried to set an indicator `@variable(model, v[1:(T+1)], Bin)` and add constraints:

```julia
@constraint(model, [t=1:(T+1)], v[t] --> {υ[t]==υ_max})
@constraint(model, [t=1:(T+1)], !v[t] --> {υ[t]==0})
@constraint(model, sum(v) ≤ v_max)

```

My code gives an error `Unable to use IndicatorToMILPBridge because element 2 in the function has a non-finite domain: 0.0 + 1.0 MOI.VariableIndex(2003)` - can I check that what I’m doing above is correct?

The complete example is below:

```Julia
using JuMP
using Ipopt
using Juniper
using Plots

## Settings

β = 0.5 # infectivity rate
γ = 0.25 # recovery rate
υ_max = 0.5 # maximum intervention
v_max = 10.0 # maximum duration of an intervention
silent = false
t0 = 0.0 # start time
tf = 100.0 # final time
δt = 0.1 # timestep
T = Int(tf/δt) # number of timesteps
S₀ = 0.99 # initial susceptible population
I₀ = 0.01 # initial infected population

## Set up model

ipopt = optimizer_with_attributes(Ipopt.Optimizer, "print_level"=>0)
optimizer = optimizer_with_attributes(Juniper.Optimizer, "nl_solver"=>ipopt)
model = Model(optimizer)

## Declare variables

@variable(model, S[1:(T+1)])
@variable(model, I[1:(T+1)])
@variable(model, υ[1:(T+1)])
@variable(model, v[1:(T+1)], Bin)

## Initial conditions
@constraint(model, S[1]==S₀)
@constraint(model, I[1]==I₀)

## Constraints on variables
@constraint(model, [t=2:(T+1)], 0 ≤ S[t] ≤ 1)
@constraint(model, [t=2:(T+1)], 0 ≤ I[t] ≤ 1)

## Constraints on control parameters
@constraint(model, [t=1:(T+1)], v[t] --> {υ[t]==υ_max})
@constraint(model, [t=1:(T+1)], !v[t] --> {υ[t]==0})
@constraint(model, δt*sum(v) ≤ v_max)

## Define auxiliary variables for infection and recovery
@NLexpression(model, infection[t=1:T], (1-exp(-(1 - υ[t]) * β * I[t] * δt)) * S[t])
@NLexpression(model, recovery[t=1:T], (1-exp(-γ*δt)) * I[t]);

## Set up timesteps
@NLconstraint(model, [t=1:T], S[t+1] == S[t] - infection[t])
@NLconstraint(model, [t=1:T], I[t+1] == I[t] + infection[t] - recovery[t])

## Minimize maximum of infection by defining an auxiliary variable
@variable(model, I_peak)
@constraint(model, 0.0 ≤ I_peak ≤ 1.0)
@constraint(model, [t=1:(T+1)], I[t] ≤ I_peak)

## Objective function
@objective(model, Min, I_peak)

if silent
    set_silent(model)
end
optimize!(model)

```

---

<div class="post-metadata">

### Author: ![blegat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/blegat/32/217090_2.png) [@blegat](https://discourse.julialang.org/u/blegat)
#### Post date: [July 11, 2024, 9:32pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/5 "2024-07-11T21:32:53Z")

</div>

The indicator constraints is reformulated using a big-M formulation so it needs bounds on the variables `u`, e.g.

```julia
@variable(model, ... <= υ[1:(T+1)] <= ...)

```

---

<div class="post-metadata">

### Author: ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)
#### Post date: [July 11, 2024, 9:48pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/6 "2024-07-11T21:48:51Z")

</div>

Thanks! How do I specify the bounds on the vector of variables? I tried `@variable(model, 0.0 ≤ υ[1:(T+1)] ≤ 1.0)` when declaring, but now I get an error `ERROR: BoundsError: attempt to access 4005-element Vector{Float64} at index [1:6007]`. How does setting the bounds when declaring a variable differ from setting a constraint later on?

---

<div class="post-metadata">

### Author: ![blegat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/blegat/32/217090_2.png) [@blegat](https://discourse.julialang.org/u/blegat)
#### Post date: [July 12, 2024, 11:56am UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/7 "2024-07-12T11:56:09Z")

</div>

Weird, it’s working for me

```julia
julia> model = Model()
A JuMP Model
Feasibility problem with:
Variables: 0
Model mode: AUTOMATIC
CachingOptimizer state: NO_OPTIMIZER
Solver name: No optimizer attached.

julia> @variable(model, 0.0 ≤ υ[1:(T+1)] ≤ 1.0)
4-element Vector{VariableRef}:
 υ[1]
 υ[2]
 υ[3]
 υ[4]

```

If you write `@constraint(model, 0 <= u[1] <= 1)` it is sent to the solver as an affine constraint, not as bounds on variables. You can use `set_lower_bound` and `set_upper_bound` to set bounds after the construction of variables.  
When rewriting Indicator constraints using Big-M formulation, we don’t compute bounds on the variables using the affine constraints, we only use the variable bounds explicitly set by the user.

---

<div class="post-metadata">

### Author: ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)
#### Post date: [July 12, 2024, 4:34pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/8 "2024-07-12T16:34:31Z")

</div>

Sorry, I wasn’t being clear - I can declare the constraints on the variable, but I now get an error. Here is the code:

```Julia
using JuMP
using Ipopt
using Juniper

β = 0.5 # infectivity rate
γ = 0.25 # recovery rate
υ_max = 0.5 # maximum intervention
v_max = 10.0 # maximum duration of an intervention
silent = false
t0 = 0.0 # start time
tf = 100.0 # final time
δt = 0.1 # timestep
T = Int(tf/δt) # number of timesteps
S₀ = 0.99 # initial susceptible population
I₀ = 0.01 # initial infected population

## Set up model

ipopt = optimizer_with_attributes(Ipopt.Optimizer, "print_level"=>0)
optimizer = optimizer_with_attributes(Juniper.Optimizer, "nl_solver"=>ipopt)
model = Model(optimizer)

## Declare variables

@variable(model, 0 ≤ S[1:(T+1)] ≤ 1)
@variable(model, 0 ≤ I[1:(T+1)] ≤ 1)
@variable(model, 0 ≤ υ[1:(T+1)] ≤ υ_max)
@variable(model, v[1:(T+1)], Bin)

## Initial conditions
@constraint(model, S[1]==S₀)
@constraint(model, I[1]==I₀)

## Constraints on control parameters
@constraint(model, [t=1:(T+1)], v[t] --> {υ[t] ≤ υ_max})
@constraint(model, [t=1:(T+1)], !v[t] --> {υ[t]==0})
@constraint(model, δt*sum(v) ≤ v_max)

## Define auxiliary variables for infection and recovery
@NLexpression(model, infection[t=1:T], (1-exp(-(1 - υ[t]) * β * I[t] * δt)) * S[t])
@NLexpression(model, recovery[t=1:T], (1-exp(-γ*δt)) * I[t]);

## Set up timesteps
@NLconstraint(model, [t=1:T], S[t+1] == S[t] - infection[t])
@NLconstraint(model, [t=1:T], I[t+1] == I[t] + infection[t] - recovery[t])

## Minimize maximum of infection by defining an auxiliary variable
@variable(model, 0 ≤ I_peak ≤ 1)
@constraint(model, [t=1:(T+1)], I[t] ≤ I_peak)

## Objective function
@objective(model, Min, I_peak)
optimize!(model)

```

The model only has 4005 parameters, but it looks like the model is trying to access more:

```Julia
nl_solver : MathOptInterface.OptimizerWithAttributes(Ipopt.Optimizer, Pair{MathOptInterface.AbstractOptimizerAttribute, Any}[MathOptInterface.RawOptimizerAttribute("print_level") => 0])
feasibility_pump : false
log_levels : [:Options, :Table, :Info]

#Variables: 4005
#IntBinVar: 1001
Obj Sense: Min

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

ERROR: BoundsError: attempt to access 4005-element Vector{Float64} at index [1:6007]
Stacktrace:
  [1] throw_boundserror(A::Vector{Float64}, I::Tuple{UnitRange{Int64}})
    @ Base ./abstractarray.jl:737
  [2] checkbounds
    @ ./abstractarray.jl:702 [inlined]
  [3] _copyto_impl!(dest::Vector{Float64}, doffs::Int64, src::Vector{Float64}, soffs::Int64, n::Int64)
    @ Base ./array.jl:374
  [4] copyto!
    @ ./array.jl:368 [inlined]
  [5] copyto!
    @ ./array.jl:388 [inlined]
  [6] _reverse_mode(d::MathOptInterface.Nonlinear.ReverseAD.NLPEvaluator, x::Vector{Float64})
    @ MathOptInterface.Nonlinear.ReverseAD ~/.julia/packages/MathOptInterface/aJZbq/src/Nonlinear/ReverseAD/reverse_mode.jl:57
  [7] eval_constraint_jacobian(d::MathOptInterface.Nonlinear.ReverseAD.NLPEvaluator, J::SubArray{…}, x::Vector{…})
    @ MathOptInterface.Nonlinear.ReverseAD ~/.julia/packages/MathOptInterface/aJZbq/src/Nonlinear/ReverseAD/mathoptinterface_api.jl:226
  [8] eval_constraint_jacobian(evaluator::MathOptInterface.Nonlinear.Evaluator{…}, J::SubArray{…}, x::Vector{…})
    @ MathOptInterface.Nonlinear ~/.julia/packages/MathOptInterface/aJZbq/src/Nonlinear/evaluator.jl:165
  [9] eval_constraint_jacobian(model::Ipopt.Optimizer, values::Vector{Float64}, x::Vector{Float64})
    @ Ipopt ~/.julia/packages/Ipopt/bqp63/src/MOI_wrapper.jl:750
 [10] (::Ipopt.var"#eval_jac_g_cb#6"{…})(x::Vector{…}, rows::Vector{…}, cols::Vector{…}, values::Vector{…})
    @ Ipopt ~/.julia/packages/Ipopt/bqp63/src/MOI_wrapper.jl:828
 [11] _Eval_Jac_G_CB(n::Int32, x_ptr::Ptr{…}, ::Int32, ::Int32, nele_jac::Int32, iRow::Ptr{…}, jCol::Ptr{…}, values_ptr::Ptr{…}, user_data::Ptr{…})
    @ Ipopt ~/.julia/packages/Ipopt/bqp63/src/C_wrapper.jl:0
 [12] IpoptSolve(prob::IpoptProblem)
    @ Ipopt ~/.julia/packages/Ipopt/bqp63/src/C_wrapper.jl:442
 [13] optimize!(model::Ipopt.Optimizer)
    @ Ipopt ~/.julia/packages/Ipopt/bqp63/src/MOI_wrapper.jl:962
 [14] optimize!(b::MathOptInterface.Bridges.LazyBridgeOptimizer{Ipopt.Optimizer})
    @ MathOptInterface.Bridges ~/.julia/packages/MathOptInterface/aJZbq/src/Bridges/bridge_optimizer.jl:367
 [15] solve_root_incumbent_model(jp::Juniper.JuniperProblem)
    @ Juniper ~/.julia/packages/Juniper/HBPrQ/src/model.jl:58
 [16] optimize!(model::Juniper.Optimizer)
    @ Juniper ~/.julia/packages/Juniper/HBPrQ/src/MOI_wrapper/MOI_wrapper.jl:301
 [17] optimize!
    @ ~/.julia/packages/MathOptInterface/aJZbq/src/Bridges/bridge_optimizer.jl:367 [inlined]
 [18] optimize!
    @ ~/.julia/packages/MathOptInterface/aJZbq/src/MathOptInterface.jl:122 [inlined]
 [19] optimize!(m::MathOptInterface.Utilities.CachingOptimizer{…})
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/cachingoptimizer.jl:321
 [20] optimize!(model::Model; ignore_optimize_hook::Bool, _differentiation_backend::MathOptInterface.Nonlinear.SparseReverseMode, kwargs::@Kwargs{})
    @ JuMP ~/.julia/packages/JuMP/7rBNn/src/optimizer_interface.jl:595
 [21] optimize!(model::Model)
    @ JuMP ~/.julia/packages/JuMP/7rBNn/src/optimizer_interface.jl:546

```

---

<div class="post-metadata">

### Author: ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [July 13, 2024, 2:24am UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/9 "2024-07-13T02:24:11Z")

</div>

I can reproduce, so I’ll take a look: [[Nonlinear.ReverseAD] BoundsError from user-code · Issue #2523 · jump-dev/MathOptInterface.jl · GitHub](https://github.com/jump-dev/MathOptInterface.jl/issues/2523)

Here’s how I would write your model though:

```Julia
using JuMP
using Ipopt
using Juniper

β = 0.5
γ = 0.25
υ_max = 0.5
v_max = 10.0
silent = false
t0 = 0.0
tf = 100.0
δt = 0.1
T = round(Int, tf / δt)
S₀ = 0.99
I₀ = 0.01

ipopt = optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)
optimizer = optimizer_with_attributes(Juniper.Optimizer, "nl_solver" => ipopt)
model = Model(optimizer)
@variables(model, begin
    0 <= S[1:(T+1)] <= 1
    0 <= I[1:(T+1)] <= 1
    0 <= υ[1:(T+1)] <= υ_max
    v[1:(T+1)], Bin
    0 <= I_peak <= 1
end)
fix(S[1], S₀)
fix(I[1], I₀)
@expressions(model, begin
    infection[t in 1:T], (1 - exp(-(1 - υ[t]) * β * I[t] * δt)) * S[t]
    recovery[t in 1:T], (1 - exp(-γ * δt)) * I[t]
end)
@constraints(model, begin
    [t in 1:T+1], υ[t] <= υ_max * v[t]
    δt * sum(v) <= v_max
    [t in 1:T], S[t+1] == S[t] - infection[t]
    [t in 1:T], I[t+1] == I[t] + infection[t] - recovery[t]
    [t in 1:T+1], I[t] <= I_peak
end)
@objective(model, Min, I_peak)
optimize!(model)

```

---

<div class="post-metadata">

### Author: ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)
#### Post date: [July 13, 2024, 3:30am UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/10 "2024-07-13T03:30:01Z")

</div>

Thanks @odow, it’s good to know the idiomatic way to do these things. Happy to know I’m not making an obvious mistake!

---

<div class="post-metadata">

### Author: ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [July 14, 2024, 5:45am UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/11 "2024-07-14T05:45:45Z")

</div>

I have identified a bug that happens when you mix the legacy nonlinear interface with a particular set of bridges (like the indicator-to-MILP bridge that is being used here).

Fix is: [[Nonlinear.ReverseAD] fix NLPBlock and final\_touch bridges by odow · Pull Request #2524 · jump-dev/MathOptInterface.jl · GitHub](https://github.com/jump-dev/MathOptInterface.jl/pull/2524)

---

<div class="post-metadata">

### Author: ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)
#### Post date: [July 14, 2024, 6:13pm UTC](https://discourse.julialang.org/t/max-of-a-vector-as-an-objective-for-jump/116911/12 "2024-07-14T18:13:10Z")

</div>

> [@sdwfrost](#):
>
> except my objective is now the maximum element in a vector (which I want to minimize)

Out of curiosity, this is a particular case of a [low order value optimization](https://link.springer.com/article/10.1007/s10898-008-9280-3) problem. And one interesting property of these problems is that, despite the fact that they are non-differentiable, the non-differentiability is “benign,” which means that you can define the objective function in a straightforward way, and use the derivative of the objective function relative to the variable for which the function assumes the maximum value. Then you can use standard derivative-based optimization methods with convergence properties similar to usual differentiable problems.
