diff --git a/docs/src/index.md b/docs/src/index.md index c1a6bf27f..d39bdefc3 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -182,13 +182,10 @@ to add the specific wrapper packages. url = {https://doi.org/10.5281/zenodo.7738525}, year = 2023} ``` - ## Reproducibility - ```@raw html
The documentation of this SciML package was built using these direct dependencies, ``` - ```@example using Pkg # hide Pkg.status() # hide diff --git a/lib/OptimizationMOI/src/nlp.jl b/lib/OptimizationMOI/src/nlp.jl index 4401ab2b9..3d0810831 100644 --- a/lib/OptimizationMOI/src/nlp.jl +++ b/lib/OptimizationMOI/src/nlp.jl @@ -317,7 +317,7 @@ end function MOI.hessian_lagrangian_structure(evaluator::MOIOptimizationNLPEvaluator) lagh = evaluator.f.lag_h !== nothing - if evaluator.f.lag_hess_prototype !== nothing + if evaluator.f.lag_hess_prototype isa SparseMatrixCSC rows, cols, _ = findnz(evaluator.f.lag_hess_prototype) return Tuple{Int, Int}[(i, j) for (i, j) in zip(rows, cols) if i <= j] end diff --git a/lib/OptimizationMetaheuristics/src/OptimizationMetaheuristics.jl b/lib/OptimizationMetaheuristics/src/OptimizationMetaheuristics.jl index d980fcd80..be4921240 100644 --- a/lib/OptimizationMetaheuristics/src/OptimizationMetaheuristics.jl +++ b/lib/OptimizationMetaheuristics/src/OptimizationMetaheuristics.jl @@ -107,9 +107,9 @@ function SciMLBase.__solve(cache::OptimizationCache{ maxiters = Optimization._check_and_convert_maxiters(cache.solver_args.maxiters) maxtime = Optimization._check_and_convert_maxtime(cache.solver_args.maxtime) - f=cache.f + f = cache.f _loss = function (θ) - if isa(f,MultiObjectiveOptimizationFunction) + if isa(f, MultiObjectiveOptimizationFunction) return cache.f(θ, cache.p) else x = cache.f(θ, cache.p) diff --git a/lib/OptimizationMetaheuristics/test/runtests.jl b/lib/OptimizationMetaheuristics/test/runtests.jl index 55d04e181..448c3d88a 100644 --- a/lib/OptimizationMetaheuristics/test/runtests.jl +++ b/lib/OptimizationMetaheuristics/test/runtests.jl @@ -53,140 +53,177 @@ Random.seed!(42) @test 10 * sol.objective < l1 # Define the benchmark functions as multi-objective problems -function sphere(x) - f1 = sum(x .^ 2) - f2 = sum((x .- 2.0) .^ 2) - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end - -function rastrigin(x) - f1 = sum(x .^ 2 .- 10 .* cos.(2 .* π .* x) .+ 10) - f2 = sum((x .- 2.0) .^ 2 .- 10 .* cos.(2 .* π .* (x .- 2.0)) .+ 10) - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end + function sphere(x) + f1 = sum(x .^ 2) + f2 = sum((x .- 2.0) .^ 2) + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end -function rosenbrock(x) - f1 = sum(100 .* (x[2:end] .- x[1:end-1] .^ 2) .^ 2 .+ (x[1:end-1] .- 1) .^ 2) - f2 = sum(100 .* ((x[2:end] .- 2.0) .- (x[1:end-1] .^ 2)) .^ 2 .+ ((x[1:end-1] .- 1.0) .^ 2)) - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end + function rastrigin(x) + f1 = sum(x .^ 2 .- 10 .* cos.(2 .* π .* x) .+ 10) + f2 = sum((x .- 2.0) .^ 2 .- 10 .* cos.(2 .* π .* (x .- 2.0)) .+ 10) + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end -function ackley(x) - f1 = -20 * exp(-0.2 * sqrt(sum(x .^ 2) / length(x))) - exp(sum(cos.(2 * π .* x)) / length(x)) + 20 + ℯ - f2 = -20 * exp(-0.2 * sqrt(sum((x .- 2.0) .^ 2) / length(x))) - exp(sum(cos.(2 * π .* (x .- 2.0))) / length(x)) + 20 + ℯ - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end + function rosenbrock(x) + f1 = sum(100 .* (x[2:end] .- x[1:(end - 1)] .^ 2) .^ 2 .+ + (x[1:(end - 1)] .- 1) .^ 2) + f2 = sum(100 .* ((x[2:end] .- 2.0) .- (x[1:(end - 1)] .^ 2)) .^ 2 .+ + ((x[1:(end - 1)] .- 1.0) .^ 2)) + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end + function ackley(x) + f1 = -20 * exp(-0.2 * sqrt(sum(x .^ 2) / length(x))) - + exp(sum(cos.(2 * π .* x)) / length(x)) + 20 + ℯ + f2 = -20 * exp(-0.2 * sqrt(sum((x .- 2.0) .^ 2) / length(x))) - + exp(sum(cos.(2 * π .* (x .- 2.0))) / length(x)) + 20 + ℯ + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end -function dtlz2(x) - g = sum((x[3:end] .- 0.5) .^ 2) - f1 = (1 + g) * cos(x[1] * π / 2) * cos(x[2] * π / 2) - f2 = (1 + g) * cos(x[1] * π / 2) * sin(x[2] * π / 2) - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end + function dtlz2(x) + g = sum((x[3:end] .- 0.5) .^ 2) + f1 = (1 + g) * cos(x[1] * π / 2) * cos(x[2] * π / 2) + f2 = (1 + g) * cos(x[1] * π / 2) * sin(x[2] * π / 2) + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end -function schaffer_n2(x) - f1 = x[1]^2 - f2 = (x[1] - 2.0)^2 - gx = [0.0] - hx = [0.0] - return [f1, f2], gx, hx -end -OBJECTIVES = Dict( - "Metaheuristics.Algorithm{NSGA2} for sphere"=> [2.1903011284699687, 3.9825426762781477], - "Metaheuristics.Algorithm{NSGA3} for sphere"=> [0.36916068436590516, 8.256797942777018], - "Metaheuristics.Algorithm{SPEA2} for sphere"=> [0.6866588142724173, 7.18284015333389], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for sphere"=> [1.6659983952552437, 4.731690734657798], - "Metaheuristics.Algorithm{MOEAD_DE} for sphere"=> [1.3118335977331483, 5.478715622895562], - "Metaheuristics.Algorithm{SMS_EMOA} for sphere"=> [0.5003293369817386, 7.837151299208113], - "Metaheuristics.Algorithm{NSGA2} for rastrigin"=> [0.0, 12.0], - "Metaheuristics.Algorithm{NSGA3} for rastrigin"=> [9.754810555001253, 11.123569741993528], - "Metaheuristics.Algorithm{SPEA2} for rastrigin"=> [0.0, 12.0], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for rastrigin"=> [2.600961284360525, 3.4282466721631755], - "Metaheuristics.Algorithm{MOEAD_DE} for rastrigin"=> [2.4963842982482607, 10.377445766099369], - "Metaheuristics.Algorithm{SMS_EMOA} for rastrigin"=> [0.0, 12.0], - "Metaheuristics.Algorithm{NSGA2} for rosenbrock"=> [17.500214034475118, 586.5039366722865], - "Metaheuristics.Algorithm{NSGA3} for rosenbrock"=> [60.58413196101549, 427.34913230512063] , - "Metaheuristics.Algorithm{SPEA2} for rosenbrock"=> [37.42314302223994, 498.8799375425481], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for rosenbrock"=> [2.600961284360525, 3.4282466721631755], - "Metaheuristics.Algorithm{MOEAD_DE} for rosenbrock"=> [12.969698120217537, 642.4135236259822], - "Metaheuristics.Algorithm{SMS_EMOA} for rosenbrock"=> [61.6898556398449, 450.62433057243777], - "Metaheuristics.Algorithm{NSGA2} for ackley"=> [2.240787163704834, 5.990002878952371], - "Metaheuristics.Algorithm{NSGA3} for ackley"=> [3.408535107623966, 5.459538604033934], - "Metaheuristics.Algorithm{SPEA2} for ackley"=> [4.440892098500626e-16, 6.593599079287213], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for ackley"=> [2.600961284360525, 3.4282466721631755], - "Metaheuristics.Algorithm{MOEAD_DE} for ackley"=> [4.440892098500626e-16, 6.593599079287213], - "Metaheuristics.Algorithm{SMS_EMOA} for ackley"=> [3.370770500897429, 5.510527199861947], - "Metaheuristics.Algorithm{NSGA2} for dtlz2"=> [0.013283104966270814, 0.010808186786590583], - "Metaheuristics.Algorithm{NSGA3} for dtlz2"=> [0.013428265441897881, 0.03589930489326534], - "Metaheuristics.Algorithm{SPEA2} for dtlz2"=> [0.019006068021099495, 0.0009905093731377751], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for dtlz2"=> [2.600961284360525, 3.4282466721631755], - "Metaheuristics.Algorithm{MOEAD_DE} for dtlz2"=> [0.027075258566241527, 0.00973958317460759], - "Metaheuristics.Algorithm{SMS_EMOA} for dtlz2"=> [0.056304481489060705, 0.026075248436234502], - "Metaheuristics.Algorithm{NSGA2} for schaffer_n2"=> [1.4034569322987955, 0.6647534264038837], - "Metaheuristics.Algorithm{NSGA3} for schaffer_n2"=> [2.7987535368174363, 0.10696329884083178], - "Metaheuristics.Algorithm{SPEA2} for schaffer_n2"=> [0.0007534237111212252, 3.8909591643988075], - "Metaheuristics.Algorithm{CCMO{NSGA2}} for schaffer_n2"=> [3.632401400816196e-17, 4.9294679997494206e-17], - "Metaheuristics.Algorithm{MOEAD_DE} for schaffer_n2"=> [2.50317097527324, 0.17460592430221922], - "Metaheuristics.Algorithm{SMS_EMOA} for schaffer_n2"=> [0.4978888767998813, 1.67543922644328], + function schaffer_n2(x) + f1 = x[1]^2 + f2 = (x[1] - 2.0)^2 + gx = [0.0] + hx = [0.0] + return [f1, f2], gx, hx + end + OBJECTIVES = Dict( + "Metaheuristics.Algorithm{NSGA2} for sphere" => [ + 2.1903011284699687, 3.9825426762781477], + "Metaheuristics.Algorithm{NSGA3} for sphere" => [ + 0.36916068436590516, 8.256797942777018], + "Metaheuristics.Algorithm{SPEA2} for sphere" => [ + 0.6866588142724173, 7.18284015333389], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for sphere" => [ + 1.6659983952552437, 4.731690734657798], + "Metaheuristics.Algorithm{MOEAD_DE} for sphere" => [ + 1.3118335977331483, 5.478715622895562], + "Metaheuristics.Algorithm{SMS_EMOA} for sphere" => [ + 0.5003293369817386, 7.837151299208113], + "Metaheuristics.Algorithm{NSGA2} for rastrigin" => [0.0, 12.0], + "Metaheuristics.Algorithm{NSGA3} for rastrigin" => [ + 9.754810555001253, 11.123569741993528], + "Metaheuristics.Algorithm{SPEA2} for rastrigin" => [0.0, 12.0], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for rastrigin" => [ + 2.600961284360525, 3.4282466721631755], + "Metaheuristics.Algorithm{MOEAD_DE} for rastrigin" => [ + 2.4963842982482607, 10.377445766099369], + "Metaheuristics.Algorithm{SMS_EMOA} for rastrigin" => [0.0, 12.0], + "Metaheuristics.Algorithm{NSGA2} for rosenbrock" => [ + 17.500214034475118, 586.5039366722865], + "Metaheuristics.Algorithm{NSGA3} for rosenbrock" => [ + 60.58413196101549, 427.34913230512063], + "Metaheuristics.Algorithm{SPEA2} for rosenbrock" => [ + 37.42314302223994, 498.8799375425481], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for rosenbrock" => [ + 2.600961284360525, 3.4282466721631755], + "Metaheuristics.Algorithm{MOEAD_DE} for rosenbrock" => [ + 12.969698120217537, 642.4135236259822], + "Metaheuristics.Algorithm{SMS_EMOA} for rosenbrock" => [ + 61.6898556398449, 450.62433057243777], + "Metaheuristics.Algorithm{NSGA2} for ackley" => [ + 2.240787163704834, 5.990002878952371], + "Metaheuristics.Algorithm{NSGA3} for ackley" => [ + 3.408535107623966, 5.459538604033934], + "Metaheuristics.Algorithm{SPEA2} for ackley" => [ + 4.440892098500626e-16, 6.593599079287213], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for ackley" => [ + 2.600961284360525, 3.4282466721631755], + "Metaheuristics.Algorithm{MOEAD_DE} for ackley" => [ + 4.440892098500626e-16, 6.593599079287213], + "Metaheuristics.Algorithm{SMS_EMOA} for ackley" => [ + 3.370770500897429, 5.510527199861947], + "Metaheuristics.Algorithm{NSGA2} for dtlz2" => [ + 0.013283104966270814, 0.010808186786590583], + "Metaheuristics.Algorithm{NSGA3} for dtlz2" => [ + 0.013428265441897881, 0.03589930489326534], + "Metaheuristics.Algorithm{SPEA2} for dtlz2" => [ + 0.019006068021099495, 0.0009905093731377751], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for dtlz2" => [ + 2.600961284360525, 3.4282466721631755], + "Metaheuristics.Algorithm{MOEAD_DE} for dtlz2" => [ + 0.027075258566241527, 0.00973958317460759], + "Metaheuristics.Algorithm{SMS_EMOA} for dtlz2" => [ + 0.056304481489060705, 0.026075248436234502], + "Metaheuristics.Algorithm{NSGA2} for schaffer_n2" => [ + 1.4034569322987955, 0.6647534264038837], + "Metaheuristics.Algorithm{NSGA3} for schaffer_n2" => [ + 2.7987535368174363, 0.10696329884083178], + "Metaheuristics.Algorithm{SPEA2} for schaffer_n2" => [ + 0.0007534237111212252, 3.8909591643988075], + "Metaheuristics.Algorithm{CCMO{NSGA2}} for schaffer_n2" => [ + 3.632401400816196e-17, 4.9294679997494206e-17], + "Metaheuristics.Algorithm{MOEAD_DE} for schaffer_n2" => [ + 2.50317097527324, 0.17460592430221922], + "Metaheuristics.Algorithm{SMS_EMOA} for schaffer_n2" => [ + 0.4978888767998813, 1.67543922644328] ) # Define the testset -@testset "Multi-Objective Optimization with Various Functions and Metaheuristics" begin - # Define the problems and their bounds - problems = [ - (sphere, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), - (rastrigin, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), - (rosenbrock, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), - (ackley, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), - (dtlz2, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), - (schaffer_n2, [0.0, 0.0, 0.0], [2.0, 0.0, 0.0]) - ] - - nobjectives = 2 - npartitions = 100 - - # Define the different algorithms - algs = [ - NSGA2(), - NSGA3(), - SPEA2(), - CCMO(NSGA2(N=100, p_m=0.001)), - MOEAD_DE(gen_ref_dirs(nobjectives, npartitions), options=Options(debug=false, iterations = 250)), - SMS_EMOA() - ] - - # Run tests for each problem and algorithm - for (prob_func, lb, ub) in problems - prob_name = string(prob_func) - for alg in algs - alg_name = string(typeof(alg)) - @testset "$alg_name on $prob_name" begin - multi_obj_fun = MultiObjectiveOptimizationFunction((x, p) -> prob_func(x)) - prob = OptimizationProblem(multi_obj_fun, lb; lb = lb, ub = ub) - if (alg_name=="Metaheuristics.Algorithm{CCMO{NSGA2}}") - sol = solve(prob, alg) - else - sol = solve(prob, alg; maxiters = 100, use_initial = true) - end - - # Tests - @test !isempty(sol.minimizer) # Check that a solution was found - - # Use sol.objective to get the objective values - key = "$alg_name for $prob_name" - value = OBJECTIVES[key] - objectives = sol.objective - @test value ≈ objectives atol=0.95 + @testset "Multi-Objective Optimization with Various Functions and Metaheuristics" begin + # Define the problems and their bounds + problems = [ + (sphere, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), + (rastrigin, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), + (rosenbrock, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), + (ackley, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), + (dtlz2, [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]), + (schaffer_n2, [0.0, 0.0, 0.0], [2.0, 0.0, 0.0]) + ] + + nobjectives = 2 + npartitions = 100 + + # Define the different algorithms + algs = [ + NSGA2(), + NSGA3(), + SPEA2(), + CCMO(NSGA2(N = 100, p_m = 0.001)), + MOEAD_DE(gen_ref_dirs(nobjectives, npartitions), + options = Options(debug = false, iterations = 250)), + SMS_EMOA() + ] + + # Run tests for each problem and algorithm + for (prob_func, lb, ub) in problems + prob_name = string(prob_func) + for alg in algs + alg_name = string(typeof(alg)) + @testset "$alg_name on $prob_name" begin + multi_obj_fun = MultiObjectiveOptimizationFunction((x, p) -> prob_func(x)) + prob = OptimizationProblem(multi_obj_fun, lb; lb = lb, ub = ub) + if (alg_name == "Metaheuristics.Algorithm{CCMO{NSGA2}}") + sol = solve(prob, alg) + else + sol = solve(prob, alg; maxiters = 100, use_initial = true) + end + + # Tests + @test !isempty(sol.minimizer) # Check that a solution was found + + # Use sol.objective to get the objective values + key = "$alg_name for $prob_name" + value = OBJECTIVES[key] + objectives = sol.objective + @test value≈objectives atol=0.95 end end end diff --git a/lib/OptimizationOptimJL/src/OptimizationOptimJL.jl b/lib/OptimizationOptimJL/src/OptimizationOptimJL.jl index 4a6a6d24e..a2df2adde 100644 --- a/lib/OptimizationOptimJL/src/OptimizationOptimJL.jl +++ b/lib/OptimizationOptimJL/src/OptimizationOptimJL.jl @@ -20,8 +20,11 @@ end SciMLBase.requiresgradient(::IPNewton) = true SciMLBase.requireshessian(::IPNewton) = true SciMLBase.requiresconsjac(::IPNewton) = true -SciMLBase.requireshessian(opt::Optim.NewtonTrustRegion) = true -SciMLBase.requireshessian(opt::Optim.Newton) = true +SciMLBase.requiresconshess(::IPNewton) = true +function SciMLBase.requireshessian(opt::Union{ + Optim.Newton, Optim.NewtonTrustRegion, Optim.KrylovTrustRegion}) + true +end SciMLBase.requiresgradient(opt::Optim.Fminbox) = true function __map_optimizer_args(cache::OptimizationCache, @@ -455,7 +458,6 @@ end using PrecompileTools PrecompileTools.@compile_workload begin - function obj_f(x, p) A = p[1] b = p[2] @@ -463,10 +465,10 @@ PrecompileTools.@compile_workload begin end function solve_nonnegative_least_squares(A, b, solver) - optf = Optimization.OptimizationFunction(obj_f, Optimization.AutoForwardDiff()) - prob = Optimization.OptimizationProblem(optf, ones(size(A, 2)), (A, b), lb=zeros(size(A, 2)), ub=Inf * ones(size(A, 2))) - x = OptimizationOptimJL.solve(prob, solver, maxiters=5000, maxtime=100) + prob = Optimization.OptimizationProblem(optf, ones(size(A, 2)), (A, b), + lb = zeros(size(A, 2)), ub = Inf * ones(size(A, 2))) + x = OptimizationOptimJL.solve(prob, solver, maxiters = 5000, maxtime = 100) return x end diff --git a/src/lbfgsb.jl b/src/lbfgsb.jl index 8a055582f..4a26b2660 100644 --- a/src/lbfgsb.jl +++ b/src/lbfgsb.jl @@ -14,12 +14,20 @@ References """ @kwdef struct LBFGS m::Int = 10 + τ = 0.5 + γ = 10.0 + λmin = -1e20 + λmax = 1e20 + μmin = 0.0 + μmax = 1e20 + ϵ = 1e-8 end SciMLBase.supports_opt_cache_interface(::LBFGS) = true SciMLBase.allowsbounds(::LBFGS) = true -# SciMLBase.requiresgradient(::LBFGS) = true +SciMLBase.requiresgradient(::LBFGS) = true SciMLBase.allowsconstraints(::LBFGS) = true +SciMLBase.requiresconsjac(::LBFGS) = true function task_message_to_string(task::Vector{UInt8}) return String(task) @@ -97,13 +105,13 @@ function SciMLBase.__solve(cache::OptimizationCache{ eq_inds = [cache.lcons[i] == cache.ucons[i] for i in eachindex(cache.lcons)] ineq_inds = (!).(eq_inds) - τ = 0.5 - γ = 10.0 - λmin = -1e20 - λmax = 1e20 - μmin = 0.0 - μmax = 1e20 - ϵ = 1e-8 + τ = cache.opt.τ + γ = cache.opt.γ + λmin = cache.opt.λmin + λmax = cache.opt.λmax + μmin = cache.opt.μmin + μmax = cache.opt.μmax + ϵ = cache.opt.ϵ λ = zeros(eltype(cache.u0), sum(eq_inds)) μ = zeros(eltype(cache.u0), sum(ineq_inds)) @@ -170,7 +178,7 @@ function SciMLBase.__solve(cache::OptimizationCache{ solver_kwargs = Base.structdiff(solver_kwargs, (; lb = nothing, ub = nothing)) for i in 1:maxiters - prev_eqcons .= cons_tmp[eq_inds] + prev_eqcons .= cons_tmp[eq_inds] .- cache.lcons[eq_inds] prevβ .= copy(β) res = optimizer(_loss, aug_grad, θ, bounds; solver_kwargs..., @@ -186,15 +194,16 @@ function SciMLBase.__solve(cache::OptimizationCache{ θ = res[2] cons_tmp .= 0.0 cache.f.cons(cons_tmp, θ) - λ = max.(min.(λmax, λ .+ ρ * cons_tmp[eq_inds]), λmin) + + λ = max.(min.(λmax, λ .+ ρ * (cons_tmp[eq_inds] .- cache.lcons[eq_inds])), λmin) β = max.(cons_tmp[ineq_inds], -1 .* μ ./ ρ) μ = min.(μmax, max.(μ .+ ρ * cons_tmp[ineq_inds], μmin)) - if max(norm(cons_tmp[eq_inds], Inf), norm(β, Inf)) > + if max(norm(cons_tmp[eq_inds] .- cache.lcons[eq_inds], Inf), norm(β, Inf)) > τ * max(norm(prev_eqcons, Inf), norm(prevβ, Inf)) ρ = γ * ρ end - if norm(cons_tmp[eq_inds], Inf) < ϵ && norm(β, Inf) < ϵ + if norm((cons_tmp[eq_inds] .- cache.lcons[eq_inds]) ./ cache.lcons[eq_inds], Inf) < ϵ && norm(β, Inf) < ϵ opt_ret = ReturnCode.Success break end diff --git a/test/ADtests.jl b/test/ADtests.jl index 029acb8e6..196c6659f 100644 --- a/test/ADtests.jl +++ b/test/ADtests.jl @@ -1,4 +1,4 @@ -using Optimization, OptimizationOptimJL, OptimizationOptimisers, Test +using Optimization, OptimizationOptimJL, OptimizationMOI, Ipopt, Test using ForwardDiff, Zygote, ReverseDiff, FiniteDiff, Tracker using ModelingToolkit, Enzyme, Random @@ -6,590 +6,92 @@ x0 = zeros(2) rosenbrock(x, p = nothing) = (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2 l1 = rosenbrock(x0) -function g!(G, x) +function g!(G, x, p = nothing) G[1] = -2.0 * (1.0 - x[1]) - 400.0 * (x[2] - x[1]^2) * x[1] G[2] = 200.0 * (x[2] - x[1]^2) end -function h!(H, x) +function h!(H, x, p = nothing) H[1, 1] = 2.0 - 400.0 * x[2] + 1200.0 * x[1]^2 H[1, 2] = -400.0 * x[1] H[2, 1] = -400.0 * x[1] H[2, 2] = 200.0 end -G1 = Array{Float64}(undef, 2) -G2 = Array{Float64}(undef, 2) -H1 = Array{Float64}(undef, 2, 2) -H2 = Array{Float64}(undef, 2, 2) +@testset "No AD" begin + optf = OptimizationFunction(rosenbrock; grad = g!, hess = h!) -g!(G1, x0) -h!(H1, x0) + prob = OptimizationProblem(optf, x0) + sol = solve(prob, Optimization.LBFGS()) -cons = (res, x, p) -> (res .= [x[1]^2 + x[2]^2]) -optf = OptimizationFunction(rosenbrock, Optimization.AutoModelingToolkit(), cons = cons) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoModelingToolkit(), - nothing, 1) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 -res = Array{Float64}(undef, 1) -optprob.cons(res, x0) -@test res == [0.0] -J = Array{Float64}(undef, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test J == [10.0, 6.0] -H3 = [Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 == [[2.0 0.0; 0.0 2.0]] + @test 10 * sol.objective < l1 + @test sol.retcode == ReturnCode.Success -function con2_c(res, x, p) - res .= [x[1]^2 + x[2]^2, x[2] * sin(x[1]) - x[1]] + sol = solve(prob, Optim.Newton()) + @test 10 * sol.objective < l1 + @test sol.retcode == ReturnCode.Success end -optf = OptimizationFunction(rosenbrock, Optimization.AutoModelingToolkit(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoModelingToolkit(), - nothing, 2) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] -G2 = Array{Float64}(undef, 2) -H2 = Array{Float64}(undef, 2, 2) +@testset "No constraint" begin + for adtype in [AutoEnzyme(), AutoForwardDiff(), AutoZygote(), AutoReverseDiff(), + AutoFiniteDiff(), AutoModelingToolkit(), AutoSparseForwardDiff(), + AutoSparseReverseDiff(), AutoSparse(AutoZygote()), AutoModelingToolkit(true, true)] + optf = OptimizationFunction(rosenbrock, adtype) -if VERSION >= v"1.9" - optf = OptimizationFunction(rosenbrock, Optimization.AutoEnzyme(), cons = cons) - optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoEnzyme(), - nothing, 1) - optprob.grad(G2, x0) - @test G1 == G2 - optprob.hess(H2, x0) - @test H1 == H2 - res = Array{Float64}(undef, 1) - optprob.cons(res, x0) - @test res == [0.0] - J = Array{Float64}(undef, 2) - optprob.cons_j(J, [5.0, 3.0]) - @test J == [10.0, 6.0] - H3 = [Array{Float64}(undef, 2, 2)] - optprob.cons_h(H3, x0) - @test H3 == [[2.0 0.0; 0.0 2.0]] + prob = OptimizationProblem(optf, x0) - G2 = Array{Float64}(undef, 2) - H2 = Array{Float64}(undef, 2, 2) + sol = solve(prob, Optim.BFGS()) + @test 10 * sol.objective < l1 + if adtype != AutoFiniteDiff() + @test sol.retcode == ReturnCode.Success + end - optf = OptimizationFunction(rosenbrock, Optimization.AutoEnzyme(), cons = con2_c) - optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoEnzyme(), - nothing, 2) - optprob.grad(G2, x0) - @test G1 == G2 - optprob.hess(H2, x0) - @test H1 == H2 - res = Array{Float64}(undef, 2) - optprob.cons(res, x0) - @test res == [0.0, 0.0] - J = Array{Float64}(undef, 2, 2) - optprob.cons_j(J, [5.0, 3.0]) - @test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) - H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] - optprob.cons_h(H3, x0) - @test H3 == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] -end - -G2 = Array{Float64}(undef, 2) -H2 = Array{Float64}(undef, 2, 2) - -optf = OptimizationFunction(rosenbrock, Optimization.AutoReverseDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoReverseDiff(), - nothing, 2) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -H3 == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -G2 = Array{Float64}(undef, 2) -H2 = Array{Float64}(undef, 2, 2) - -optf = OptimizationFunction(rosenbrock, Optimization.AutoReverseDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, - Optimization.AutoReverseDiff(compile = true), - nothing, 2) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -H3 == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -G2 = Array{Float64}(undef, 2) -H2 = Array{Float64}(undef, 2, 2) - -optf = OptimizationFunction(rosenbrock, Optimization.AutoZygote(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoZygote(), - nothing, 2) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -H3 == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -optf = OptimizationFunction(rosenbrock, Optimization.AutoModelingToolkit(true, true), - cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, - Optimization.AutoModelingToolkit(true, true), - nothing, 2) -using SparseArrays -sH = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -@test findnz(sH)[1:2] == findnz(optprob.hess_prototype)[1:2] -optprob.hess(sH, x0) -@test sH == H2 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -sJ = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -@test findnz(sJ)[1:2] == findnz(optprob.cons_jac_prototype)[1:2] -optprob.cons_j(sJ, [5.0, 3.0]) -@test all(isapprox(sJ, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -sH3 = [sparse([1, 2], [1, 2], zeros(2)), sparse([1, 1, 2], [1, 2, 1], zeros(3))] -@test getindex.(findnz.(sH3), Ref([1, 2])) == - getindex.(findnz.(optprob.cons_hess_prototype), Ref([1, 2])) -optprob.cons_h(sH3, x0) -@test Array.(sH3) == [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -optf = OptimizationFunction(rosenbrock, Optimization.AutoForwardDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoForwardDiff(), - nothing) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 ≈ H2 - -prob = OptimizationProblem(optf, x0) - -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 -@test sol.retcode == ReturnCode.Success - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 -@test sol.retcode == ReturnCode.Success - -sol = solve(prob, Optim.KrylovTrustRegion()) -@test 10 * sol.objective < l1 -@test sol.retcode == ReturnCode.Success - -optf = OptimizationFunction(rosenbrock, Optimization.AutoZygote()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoZygote(), nothing) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 - -prob = OptimizationProblem(optf, x0) - -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.KrylovTrustRegion()) -@test 10 * sol.objective < l1 - -optf = OptimizationFunction(rosenbrock, Optimization.AutoReverseDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoReverseDiff(), - nothing) -optprob.grad(G2, x0) -@test G1 == G2 -optprob.hess(H2, x0) -@test H1 == H2 - -prob = OptimizationProblem(optf, x0) -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.KrylovTrustRegion()) -@test 10 * sol.objective < l1 - -optf = OptimizationFunction(rosenbrock, Optimization.AutoTracker()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoTracker(), nothing) -optprob.grad(G2, x0) -@test G1 == G2 -@test_broken optprob.hess(H2, x0) - -prob = OptimizationProblem(optf, x0) - -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -@test_broken solve(prob, Newton()) - -optf = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoFiniteDiff(), - nothing) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-6 - -prob = OptimizationProblem(optf, x0) -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 + sol = solve(prob, Optim.Newton()) + @test 10 * sol.objective < l1 + @test sol.retcode == ReturnCode.Success -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 + sol = solve(prob, Optim.KrylovTrustRegion()) + @test 10 * sol.objective < l1 + @test sol.retcode == ReturnCode.Success -sol = solve(prob, Optim.KrylovTrustRegion()) -@test sol.objective < l1 #the loss doesn't go below 5e-1 here - -sol = solve(prob, Optimisers.Adam(0.1), maxiters = 1000) -@test 10 * sol.objective < l1 - -# Test new constraints -cons = (res, x, p) -> (res .= [x[1]^2 + x[2]^2]) -optf = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(), cons = cons) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoFiniteDiff(), - nothing, 1) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-6 -res = Array{Float64}(undef, 1) -optprob.cons(res, x0) -@test res == [0.0] -optprob.cons(res, [1.0, 4.0]) -@test res == [17.0] -J = zeros(1, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test J ≈ [10.0 6.0] -H3 = [Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0]] - -# H4 = Array{Float64}(undef, 2, 2) -# μ = randn(1) -# σ = rand() -# optprob.lag_h(H4, x0, σ, μ) -# @test H4≈σ * H1 + μ[1] * H3[1] rtol=1e-6 - -cons_jac_proto = Float64.(sparse([1 1])) # Things break if you only use [1 1]; see FiniteDiff.jl -cons_jac_colors = 1:2 -optf = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(), cons = cons, - cons_jac_prototype = cons_jac_proto, - cons_jac_colorvec = cons_jac_colors) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoFiniteDiff(), - nothing, 1) -@test optprob.cons_jac_prototype == sparse([1.0 1.0]) # make sure it's still using it -@test optprob.cons_jac_colorvec == 1:2 -J = zeros(1, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test J ≈ [10.0 6.0] - -function con2_c(res, x, p) - res .= [x[1]^2 + x[2]^2, x[2] * sin(x[1]) - x[1]] + sol = solve(prob, Optimization.LBFGS(), maxiters = 1000) + @test 10 * sol.objective < l1 + @test sol.retcode == ReturnCode.Success + end end -optf = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoFiniteDiff(), - nothing, 2) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-6 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res == [0.0, 0.0] -optprob.cons(res, [1.0, 2.0]) -@test res ≈ [5.0, 0.682941969615793] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -cons_jac_proto = Float64.(sparse([1 1; 1 1])) -cons_jac_colors = 1:2 -optf = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(), cons = con2_c, - cons_jac_prototype = cons_jac_proto, - cons_jac_colorvec = cons_jac_colors) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoFiniteDiff(), - nothing, 2) -@test optprob.cons_jac_prototype == sparse([1.0 1.0; 1.0 1.0]) # make sure it's still using it -@test optprob.cons_jac_colorvec == 1:2 -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H2 = Array{Float64}(undef, 2, 2) -optprob.hess(H2, [5.0, 3.0]) -@test all(isapprox(H2, [28802.0 -2000.0; -2000.0 200.0]; rtol = 1e-3)) -cons_j = (J, θ, p) -> optprob.cons_j(J, θ) -hess = (H, θ, p) -> optprob.hess(H, θ) -sH = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -sJ = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -optf = OptimizationFunction(rosenbrock, Optimization.AutoForwardDiff(), hess = hess, - hess_prototype = copy(sH), cons = con2_c, cons_j = cons_j, - cons_jac_prototype = copy(sJ)) -optprob1 = Optimization.instantiate_function(optf, x0, Optimization.AutoForwardDiff(), - nothing, 2) -@test optprob1.hess_prototype == sparse([0.0 0.0; 0.0 0.0]) # make sure it's still using it -optprob1.hess(sH, [5.0, 3.0]) -@test all(isapprox(sH, [28802.0 -2000.0; -2000.0 200.0]; rtol = 1e-3)) -@test optprob1.cons_jac_prototype == sparse([0.0 0.0; 0.0 0.0]) # make sure it's still using it -optprob1.cons_j(sJ, [5.0, 3.0]) -@test all(isapprox(sJ, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) +@testset "One constraint" begin + for adtype in [AutoEnzyme(), AutoForwardDiff(), AutoZygote(), AutoReverseDiff(), AutoModelingToolkit(), AutoSparseForwardDiff(), AutoSparseReverseDiff(), AutoSparse(AutoZygote()), AutoModelingToolkit(true, true)] + cons = (res, x, p) -> (res[1] = x[1]^2 + x[2]^2 - 1.0; return nothing) + optf = OptimizationFunction(rosenbrock, adtype, cons = cons) + @show adtype + prob = OptimizationProblem(optf, x0, lb = [-1.0, -1.0], ub = [1.0, 1.0], lcons = [0.0], ucons = [0.0]) -grad = (G, θ, p) -> optprob.grad(G, θ) -hess = (H, θ, p) -> optprob.hess(H, θ) -cons_j = (J, θ, p) -> optprob.cons_j(J, θ) -cons_h = (res, θ, p) -> optprob.cons_h(res, θ) -sH = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -sJ = sparse([1, 1, 2, 2], [1, 2, 1, 2], zeros(4)) -sH3 = [sparse([1, 2], [1, 2], zeros(2)), sparse([1, 1, 2], [1, 2, 1], zeros(3))] -optf = OptimizationFunction(rosenbrock, SciMLBase.NoAD(), grad = grad, hess = hess, - cons = con2_c, cons_j = cons_j, cons_h = cons_h, - hess_prototype = sH, cons_jac_prototype = sJ, - cons_hess_prototype = sH3) -optprob2 = Optimization.instantiate_function(optf, x0, SciMLBase.NoAD(), nothing, 2) -optprob2.hess(sH, [5.0, 3.0]) -@test all(isapprox(sH, [28802.0 -2000.0; -2000.0 200.0]; rtol = 1e-3)) -optprob2.cons_j(sJ, [5.0, 3.0]) -@test all(isapprox(sJ, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -optprob2.cons_h(sH3, [5.0, 3.0]) -@test sH3 ≈ [ - [2.0 0.0; 0.0 2.0], - [2.8767727327346804 0.2836621681849162; 0.2836621681849162 -6.622738308376736e-9] -] + sol = solve(prob, Optimization.LBFGS(), maxiters = 1000) + @test 10 * sol.objective < l1 -# Can we solve problems? Using AutoForwardDiff to test since we know that works -for consf in [cons, con2_c] - optf1 = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(); cons = consf) - lcons = consf == cons ? [0.2] : [0.2, -0.81] - ucons = consf == cons ? [0.55] : [0.55, -0.1] - prob1 = OptimizationProblem(optf1, [0.3, 0.5], lb = [0.2, 0.4], ub = [0.6, 0.8], - lcons = lcons, ucons = ucons) - sol1 = solve(prob1, Optim.IPNewton()) - @test sol1.retcode == ReturnCode.Success - optf2 = OptimizationFunction(rosenbrock, Optimization.AutoForwardDiff(); cons = consf) - prob2 = OptimizationProblem(optf2, [0.3, 0.5], lb = [0.2, 0.4], ub = [0.6, 0.8], - lcons = lcons, ucons = ucons) - sol2 = solve(prob2, Optim.IPNewton()) - @test sol2.retcode == ReturnCode.Success - @test sol1.objective≈sol2.objective rtol=1e-4 - @test sol1.u ≈ sol2.u - res = Array{Float64}(undef, length(lcons)) - consf(res, sol1.u, nothing) - @test lcons[1] ≤ res[1] ≤ ucons[1] - if consf == con2_c - @test lcons[2] ≤ res[2] ≤ ucons[2] - end - - lcons = consf == cons ? [0.2] : [0.2, 0.5] - ucons = consf == cons ? [0.2] : [0.2, 0.5] - optf1 = OptimizationFunction(rosenbrock, Optimization.AutoFiniteDiff(); cons = consf) - prob1 = OptimizationProblem(optf1, [0.5, 0.5], lcons = lcons, ucons = ucons) - sol1 = solve(prob1, Optim.IPNewton()) - @test sol1.retcode == ReturnCode.Success - optf2 = OptimizationFunction(rosenbrock, Optimization.AutoForwardDiff(); cons = consf) - prob2 = OptimizationProblem(optf2, [0.5, 0.5], lcons = lcons, ucons = ucons) - sol2 = solve(prob2, Optim.IPNewton()) - @test sol2.retcode == ReturnCode.Success - @test sol1.objective≈sol2.objective rtol=1e-4 - @test sol1.u≈sol2.u rtol=1e-4 - res = Array{Float64}(undef, length(lcons)) - consf(res, sol1.u, nothing) - @test res[1]≈lcons[1] rtol=1e-1 - if consf == con2_c - @test res[2]≈lcons[2] rtol=1e-2 + sol = solve(prob, Ipopt.Optimizer(), max_iter = 1000; print_level = 0) + @test 10 * sol.objective < l1 end end -using SparseDiffTools - -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseFiniteDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseFiniteDiff(), - nothing, 2) -G2 = Array{Float64}(undef, 2) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-4 -H2 = Array{Float64}(undef, 2, 2) -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-4 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res≈[0.0, 0.0] atol=1e-4 -optprob.cons(res, [1.0, 2.0]) -@test res ≈ [5.0, 0.682941969615793] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseFiniteDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseFiniteDiff(), - nothing) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-4 - -prob = OptimizationProblem(optf, x0) -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 - -Random.seed!(1234) -#at 0,0 it gives error because of the inaccuracy of the hessian and hv calculations -prob = OptimizationProblem(optf, x0 + rand(2)) -sol = solve(prob, Optim.KrylovTrustRegion()) -@test sol.objective < l1 +@testset "Two constraints" begin + for adtype in [AutoEnzyme(), AutoForwardDiff(), AutoZygote(), AutoReverseDiff(), + AutoModelingToolkit(), AutoSparseForwardDiff(), + AutoSparseReverseDiff(), AutoSparse(AutoZygote()), AutoModelingToolkit(true, true)] + function con2_c(res, x, p) + res[1] = x[1]^2 + x[2]^2 + res[2] = x[2] * sin(x[1]) - x[1] + return nothing + end + optf = OptimizationFunction(rosenbrock, adtype, cons = con2_c) -sol = solve(prob, Optimisers.ADAM(0.1), maxiters = 1000) -@test 10 * sol.objective < l1 + prob = OptimizationProblem(optf, x0, lb = [-1.0, -1.0], ub = [1.0, 1.0], lcons = [1.0, -2.0], ucons = [1.0, 2.0]) -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseForwardDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseForwardDiff(), - nothing, 2) -G2 = Array{Float64}(undef, 2) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-4 -H2 = Array{Float64}(undef, 2, 2) -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-4 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res≈[0.0, 0.0] atol=1e-4 -optprob.cons(res, [1.0, 2.0]) -@test res ≈ [5.0, 0.682941969615793] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] + sol = solve(prob, Optimization.LBFGS(), maxiters = 1000) + @test 10 * sol.objective < l1 -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseForwardDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseForwardDiff(), - nothing) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-6 - -prob = OptimizationProblem(optf, x0) -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.KrylovTrustRegion()) -@test sol.objective < l1 - -sol = solve(prob, Optimisers.ADAM(0.1), maxiters = 1000) -@test 10 * sol.objective < l1 - -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseReverseDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, - Optimization.AutoSparseReverseDiff(true), - nothing, 2) -G2 = Array{Float64}(undef, 2) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-4 -H2 = Array{Float64}(undef, 2, 2) -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-4 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res≈[0.0, 0.0] atol=1e-4 -optprob.cons(res, [1.0, 2.0]) -@test res ≈ [5.0, 0.682941969615793] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseReverseDiff(), cons = con2_c) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseReverseDiff(), - nothing, 2) -G2 = Array{Float64}(undef, 2) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-4 -H2 = Array{Float64}(undef, 2, 2) -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-4 -res = Array{Float64}(undef, 2) -optprob.cons(res, x0) -@test res≈[0.0, 0.0] atol=1e-4 -optprob.cons(res, [1.0, 2.0]) -@test res ≈ [5.0, 0.682941969615793] -J = Array{Float64}(undef, 2, 2) -optprob.cons_j(J, [5.0, 3.0]) -@test all(isapprox(J, [10.0 6.0; -0.149013 -0.958924]; rtol = 1e-3)) -H3 = [Array{Float64}(undef, 2, 2), Array{Float64}(undef, 2, 2)] -optprob.cons_h(H3, x0) -@test H3 ≈ [[2.0 0.0; 0.0 2.0], [-0.0 1.0; 1.0 0.0]] - -optf = OptimizationFunction(rosenbrock, Optimization.AutoSparseReverseDiff()) -optprob = Optimization.instantiate_function(optf, x0, Optimization.AutoSparseReverseDiff(), - nothing) -optprob.grad(G2, x0) -@test G1≈G2 rtol=1e-6 -optprob.hess(H2, x0) -@test H1≈H2 rtol=1e-6 - -prob = OptimizationProblem(optf, x0) -sol = solve(prob, Optim.BFGS()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.Newton()) -@test 10 * sol.objective < l1 - -sol = solve(prob, Optim.KrylovTrustRegion()) -@test sol.objective < l1 - -sol = solve(prob, Optimisers.ADAM(0.1), maxiters = 1000) -@test 10 * sol.objective < l1 + sol = solve(prob, Ipopt.Optimizer(), max_iter = 1000; print_level = 0) + @test 10 * sol.objective < l1 + end +end diff --git a/test/Project.toml b/test/Project.toml index af64fbe0c..ba5236ded 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -8,6 +8,7 @@ Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" FiniteDiff = "6a86dc24-6348-571c-b903-95158fe2bd41" Flux = "587475ba-b771-5e3f-ad9e-33799f191a9c" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" +Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9" IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e" Lux = "b2108857-7c20-44ae-9111-449ecde12c47" ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" diff --git a/test/runtests.jl b/test/runtests.jl index 8b85ae3a5..0bf8d63f2 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -14,8 +14,11 @@ function activate_subpkg_env(subpkg) Pkg.instantiate() end -if GROUP == "All" || GROUP == "Core" || GROUP == "GPU" || - GROUP == "OptimizationPolyalgorithms" +if GROUP == "All" || GROUP == "Core" + dev_subpkg("OptimizationOptimJL") + dev_subpkg("OptimizationOptimisers") + dev_subpkg("OptimizationMOI") +elseif GROUP == "GPU" || GROUP == "OptimizationPolyalgorithms" dev_subpkg("OptimizationOptimJL") dev_subpkg("OptimizationOptimisers") end diff --git a/test/stdout.txt b/test/stdout.txt new file mode 100644 index 000000000..8a263fca6 --- /dev/null +++ b/test/stdout.txt @@ -0,0 +1 @@ +ErrorException("type Array has no field nzval") \ No newline at end of file