@ChrisRackauckas can you specify what you mean precisely? In the MWE I wrote up the call I made to PRIMA (through Optimization.jl). I’ve included a version with the cobyla function from PRIMA.jl directly, which does give the right solution:
using Random, Optimization, OptimizationMOI, OptimizationNLopt, OptimizationOptimJL, OptimizationPRIMA, DifferentiationInterface, Ipopt, PRIMA
Random.seed!(42)
function RunTest()
rosenbrock(x, p) = (p[1] - x[1])^2 + p[2] * (x[2] - x[1]^2)^2
x0 = [0.5, 0.5]
_p = [1.0, 100.0]
f(x) = [x[1]^2 + x[2]^2, x[1] * x[2]]
cons(res, x, p) = (res .= f(x))
optprob = OptimizationFunction(rosenbrock, SecondOrder(AutoForwardDiff(), AutoForwardDiff()), cons = cons)
#prob = OptimizationProblem(optprob, x0, _p, lb = [0., 0.], ub = [1., 1.], lcons = [-1., -1.0], ucons = [2., 2.0])
prob = OptimizationProblem(optprob, x0, _p, lb = [0., 0.], ub = [1., 1.], lcons = [-.8, -2.0], ucons = [0.8, 2.0])
opt_vec = [COBYLA(), IPNewton()]
opt_names = ["COBYLA", "IPNewton"]
for (opt, opt_name) in zip(opt_vec, opt_names)
sol = solve(prob, opt)
print("====\n$opt_name\n$(sol)")
end
x, info = cobyla(x -> rosenbrock(x, _p), x0; nonlinear_ineq = x -> f(x) .- [0.8, 2.])
print("====\nusing cobyla directly: $x\n$info\n")
end
RunTest()
If I had to guess it may be that PRIMA.jl requires constraints to be of the “f <= 0” type rather than of the “lc <= f <= uc” type, perhaps there is some error in translation here when Optimization.jl wraps COBYLA() and calls cobyla() in PRIMA.jl? I’m not too much aware of the exact implementation here.
====
COBYLA
retcode: Success
u: [1.58113883008419e-6, 4.743416490252569e-6]
Final objective value: 0.9999968399748375
====
IPNewton
retcode: Success
u: [0.7247018404470709, 0.5242206047577943]
Final objective value: 0.07588358475785834
====
using cobyla directly: [0.7247016181219436, 0.5242209192620982]
PRIMA.Info(0.07588358340753011, 365, PRIMA.SMALL_TR_RADIUS, 7.500562615447848e-9, Float64[], [7.500562615447848e-9, -1.6200962515573847])