Max of a vector as an objective for JuMP?

Specific Domains Optimization (Mathematical)
8

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

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:

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
1 Like