How to add a linear constraint to optimize() using optim.jl?

Specific Domains Optimization (Mathematical)
1

Hi,

I wanted to add a linear constraint to a maximization problem using optim.jl. It is a linear constraint and cannot be done by box constrain. As for algorithms, I will use both gradient free and Gradient required methods. Attached is a MWE.

g_guess = [0.01,0.01,0.01]
    
    #lower = [0,0,0]
    #upper = [1,1,1]
    
    #func =  TwiceDifferentiable(g -> -ll_wrapper_fB_fG_p_s(g,rdata_prepro,theta,options),
    #                            g_guess);
    a=[0.3,-0.2,0.11]
    b=[-0.5,0.8,0.1]
    func(g) = -sum(log.(cdf.(Normal(),((1-g[1]) .* a .+ (1-g[1]-g[2]) .* b)./g[3])))
    
    opt=optimize(func,g_guess,Optim.Options(x_tol = 1e-5,f_tol =  1e-5,g_tol=1e-5, store_trace = true, show_trace=true,extended_trace=true,show_every=1))

    #estimates for theta parameters
    g_hat = Optim.minimizer(opt)
    
    print("estimates for theta are: $g_hat")
    
    print("opt info is: $opt")

I wanted to add a constraint: g[2]+g[3]<1. Is it possible in optim.jl? If not, is there any alternative ways to do that?

Thank you in advance!!

1 Like
2

Use JuMP.jl

2 Likes
3

As @jling suggests, you can use JuMP for this:

julia> using JuMP

julia> import Ipopt

julia> import Distributions

julia> function main(; a, b, g_guess)
           N = length(a)
           model = Model(Ipopt.Optimizer)
           cdf(x) = Distributions.cdf(Distributions.Normal(), x)
           register(model, :cdf, 1, cdf; autodiff = true)
           @variable(model, g[i = 1:3], start = g_guess[i])
           @NLexpression(
               model, 
               expr[i=1:N], 
               ((1 - g[1]) * a[i] + (1 - g[1] - g[2]) * b[i]) / g[3],
           )
           @NLobjective(model, Max, sum(log(cdf(expr[i])) for i in 1:N))
           @constraint(model, g[2] + g[3] <= 1)
           optimize!(model)
           print(solution_summary(model))
           return value.(g)
       end
main (generic function with 2 methods)

julia> main(a = [0.3, -0.2, 0.11], b = [-0.5, 0.8, 0.1], g_guess = [0.01, 0.01, 0.01])
This is Ipopt version 3.14.4, running with linear solver MUMPS 5.4.1.

Number of nonzeros in equality constraint Jacobian...:        0
Number of nonzeros in inequality constraint Jacobian.:        2
Number of nonzeros in Lagrangian Hessian.............:        6

Total number of variables............................:        3
                     variables with only lower bounds:        0
                variables with lower and upper bounds:        0
                     variables with only upper bounds:        0
Total number of equality constraints.................:        0
Total number of inequality constraints...............:        1
        inequality constraints with only lower bounds:        0
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:        1

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0 -1.9012671e+02 0.00e+00 6.55e+01  -1.0 0.00e+00    -  0.00e+00 0.00e+00   0
   1 -1.0776751e+02 0.00e+00 4.17e+01  -1.0 5.53e-03   2.0 1.00e+00 1.00e+00f  1
   2 -6.1174561e+01 0.00e+00 1.75e+01  -1.0 8.22e-03   1.5 1.00e+00 1.00e+00f  1
   3 -3.4854534e+01 0.00e+00 7.32e+00  -1.0 1.13e-02   1.0 1.00e+00 1.00e+00f  1
   4 -2.0115941e+01 0.00e+00 3.05e+00  -1.0 1.86e-02   0.6 1.00e+00 1.00e+00f  1
   5 -1.2175019e+01 0.00e+00 1.28e+00  -1.0 3.25e-02   0.1 1.00e+00 1.00e+00f  1
   6 -8.5215155e+00 0.00e+00 5.69e-01  -1.0 6.55e-02  -0.4 1.00e+00 1.00e+00f  1
   7 -4.3474423e+00 0.00e+00 2.04e-01  -1.7 1.29e-01  -0.9 1.00e+00 1.00e+00f  1
   8 -1.0334178e-05 0.00e+00 5.16e-02  -1.7 2.28e+00  -1.3 1.00e+00 1.27e-01f  3
   9 -1.1878174e-08 0.00e+00 5.76e-04  -2.5 3.00e-02  -1.8 1.00e+00 1.00e+00f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  10 -3.6985459e-09 0.00e+00 4.34e-08  -5.7 2.60e-03    -  1.00e+00 1.00e+00f  1
  11 -1.6667980e-08 0.00e+00 4.10e-08  -8.6 2.29e+00    -  8.84e-01 1.00e+00f  1
  12 -1.6581511e-08 0.00e+00 2.77e-08  -8.6 5.45e-06  -2.3 1.00e+00 1.00e+00h  1
  13 -2.3329098e-08 0.00e+00 1.90e-08  -8.6 3.56e+00    -  9.68e-01 1.00e+00f  1
  14 -2.3272493e-08 0.00e+00 1.30e-08  -9.0 7.68e-06  -2.8 1.00e+00 1.00e+00f  1
  15 -2.7242802e-08 0.00e+00 9.33e-09  -9.0 7.23e+00    -  9.68e-01 1.00e+00f  1

Number of Iterations....: 15

                                   (scaled)                 (unscaled)
Objective...............:   7.2942216959846420e-11   -2.7242802121052110e-08
Dual infeasibility......:   9.3335199888424968e-09    3.4859269260885378e-06
Constraint violation....:   0.0000000000000000e+00    0.0000000000000000e+00
Variable bound violation:   0.0000000000000000e+00    0.0000000000000000e+00
Complementarity.........:   1.6800458301956473e-18    6.2747141523699448e-16
Overall NLP error.......:   9.3335199888424968e-09    3.4859269260885378e-06


Number of objective function evaluations             = 19
Number of objective gradient evaluations             = 16
Number of equality constraint evaluations            = 0
Number of inequality constraint evaluations          = 19
Number of equality constraint Jacobian evaluations   = 0
Number of inequality constraint Jacobian evaluations = 1
Number of Lagrangian Hessian evaluations             = 15
Total seconds in IPOPT                               = 0.115

EXIT: Optimal Solution Found.
* Solver : Ipopt

* Status
  Termination status : LOCALLY_SOLVED
  Primal status      : FEASIBLE_POINT
  Dual status        : FEASIBLE_POINT
  Message from the solver:
  "Solve_Succeeded"

* Candidate solution
  Objective value      : -2.72428e-08

* Work counters
  Solve time (sec)   : 1.15788e-01
3-element Vector{Float64}:
 13.642843753497813
 -7.692989230060321
 -0.24194773252752094
3 Likes