Create ode.md

docs/src/optimization_packages/ode.md

@@ -0,0 +1,63 @@
+# OptimizationODE.jl
+
+**OptimizationODE.jl** provides ODE-based optimization methods as a solver plugin for [SciML's Optimization.jl](https://github.com/SciML/Optimization.jl). It wraps various ODE solvers to perform gradient-based optimization using continuous-time dynamics.
+
+## Installation
+
+```julia
+using Pkg
+Pkg.add(url="OptimizationODE.jl")
+```
+
+## Usage
+
+```julia
+using OptimizationODE, Optimization, ADTypes, SciMLBase
+
+function f(x, p)
+    return sum(abs2, x)
+end
+
+function g!(g, x, p)
+    @. g = 2 * x
+end
+
+x0 = [2.0, -3.0]
+p = []
+
+f_manual = OptimizationFunction(f, SciMLBase.NoAD(); grad = g!)
+prob_manual = OptimizationProblem(f_manual, x0)
+
+opt = ODEGradientDescent()
+sol = solve(prob_manual, opt; dt=0.01, maxiters=50_000)
+
+@show sol.u
+@show sol.objective
+```
+
+## Local Gradient-based Optimizers
+
+All provided optimizers are **gradient-based local optimizers** that solve optimization problems by integrating gradient-based ODEs to convergence:
+
+* `ODEGradientDescent()` — performs basic gradient descent using the explicit Euler method. This is a simple and efficient method suitable for small-scale or well-conditioned problems.
+
+* `RKChebyshevDescent()` — uses the ROCK2 solver, a stabilized explicit Runge-Kutta method suitable for stiff problems. It allows larger step sizes while maintaining stability.
+
+* `RKAccelerated()` — leverages the Tsit5 method, a 5th-order Runge-Kutta solver that achieves faster convergence for smooth problems by improving integration accuracy.
+
+* `HighOrderDescent()` — applies Vern7, a high-order (7th-order) explicit Runge-Kutta method for even more accurate integration. This can be beneficial for problems requiring high precision.
+
+## DAE-based Optimizers
+
+In addition to ODE-based optimizers, OptimizationODE.jl provides optimizers for differential-algebraic equation (DAE) constrained problems:
+
+* `DAEMassMatrix()` — uses the Rodas5 solver (from OrdinaryDiffEq.jl) for DAE problems with a mass matrix formulation.
+
+* `DAEIndexing()` — uses the IDA solver (from Sundials.jl) for DAE problems with index variable support.
+
+You can also define a custom optimizer using the generic `ODEOptimizer(solver)` or `DAEOptimizer(solver)` constructor by supplying any ODE or DAE solver supported by [OrdinaryDiffEq.jl](https://docs.sciml.ai/DiffEqDocs/stable/solvers/ode_solve/) or [Sundials.jl](https://github.com/SciML/Sundials.jl).
+
+## Interface Details
+
+All optimizers require gradient information (either via automatic differentiation or manually provided `grad!`). The optimization is performed by integrating the ODE defined by the negative gradient until a steady state is reached.
+

lib/OptimizationBBO/src/OptimizationBBO.jl

@@ -53,9 +53,25 @@ function __map_optimizer_args(prob::Optimization.OptimizationCache, opt::BBO;
     if !isnothing(reltol)
         @warn "common reltol is currently not used by $(opt)"
     end
-    mapped_args = (; kwargs...)
-    mapped_args = (; mapped_args..., Method = opt.method,
-        SearchRange = [(prob.lb[i], prob.ub[i]) for i in 1:length(prob.lb)])
+
+    # Determine number of objectives for multi-objective problems
+    if isa(prob.f, MultiObjectiveOptimizationFunction)
+        num_objectives = length(prob.f.cost_prototype)
+        mapped_args = (; kwargs...)
+        mapped_args = (; mapped_args..., Method = opt.method,
+            SearchRange = [(prob.lb[i], prob.ub[i]) for i in 1:length(prob.lb)],
+            NumDimensions = length(prob.lb),
+            NumObjectives = num_objectives)
+        # FitnessScheme should be in opt, not the function
+        if hasproperty(opt, :FitnessScheme)
+            mapped_args = (; mapped_args..., FitnessScheme = opt.FitnessScheme)
+        end
+    else
+        mapped_args = (; kwargs...)
+        mapped_args = (; mapped_args..., Method = opt.method,
+            SearchRange = [(prob.lb[i], prob.ub[i]) for i in 1:length(prob.lb)],
+            NumDimensions = length(prob.lb))
+    end
 
     if !isnothing(callback)
         mapped_args = (; mapped_args..., CallbackFunction = callback,
@@ -144,8 +160,16 @@ function SciMLBase.__solve(cache::Optimization.OptimizationCache{
     maxiters = Optimization._check_and_convert_maxiters(cache.solver_args.maxiters)
     maxtime = Optimization._check_and_convert_maxtime(cache.solver_args.maxtime)
 
-    _loss = function (θ)
-        cache.f(θ, cache.p)
+
+    # Multi-objective: use out-of-place or in-place as appropriate
+    if isa(cache.f, MultiObjectiveOptimizationFunction)
+        if cache.f.iip
+            _loss = θ -> (cost = similar(cache.f.cost_prototype); cache.f.f(cost, θ, cache.p); cost)
+        else
+            _loss = θ -> cache.f.f(θ, cache.p)
+        end
+    else
+        _loss = θ -> cache.f(θ, cache.p)
     end
 
     opt_args = __map_optimizer_args(cache, cache.opt;

lib/OptimizationBBO/test/runtests.jl

@@ -34,6 +34,50 @@ using Test
         push!(loss_history, fitness)
         return false
     end
+
+    @testset "In-place Multi-Objective Optimization" begin
+        function inplace_multi_obj!(cost, x, p)
+            cost[1] = sum(x .^ 2)
+            cost[2] = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return nothing
+        end
+        u0 = [0.25, 0.25]
+        lb = [0.0, 0.0]
+        ub = [2.0, 2.0]
+        cost_prototype = zeros(2)
+        mof_inplace = MultiObjectiveOptimizationFunction(inplace_multi_obj!; cost_prototype=cost_prototype)
+        prob_inplace = Optimization.OptimizationProblem(mof_inplace, u0; lb=lb, ub=ub)
+        sol_inplace = solve(prob_inplace, opt, NumDimensions=2, FitnessScheme=ParetoFitnessScheme{2}(is_minimizing=true))
+        @test sol_inplace ≠ nothing
+        @test length(sol_inplace.objective) == 2
+        @test sol_inplace.objective[1] ≈ 6.9905986e-18 atol=1e-3
+        @test sol_inplace.objective[2] ≈ 1.7763568e-15 atol=1e-3
+    end
+
+    @testset "Custom coalesce for Multi-Objective" begin
+        function multi_obj_tuple(x, p)
+            f1 = sum(x .^ 2)
+            f2 = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return (f1, f2)
+        end
+        coalesce_sum(cost, x, p) = sum(cost)
+        mof_coalesce = MultiObjectiveOptimizationFunction(multi_obj_tuple; coalesce=coalesce_sum)
+        prob_coalesce = Optimization.OptimizationProblem(mof_coalesce, u0; lb=lb, ub=ub)
+        sol_coalesce = solve(prob_coalesce, opt, NumDimensions=2, FitnessScheme=ParetoFitnessScheme{2}(is_minimizing=true))
+        @test sol_coalesce ≠ nothing
+        @test sol_coalesce.objective[1] ≈ 6.9905986e-18 atol=1e-3
+        @test sol_coalesce.objective[2] ≈ 1.7763568e-15 atol=1e-3
+        @test mof_coalesce.coalesce([1.0, 2.0], [0.0, 0.0], nothing) == 3.0
+    end
+
+    @testset "Error if in-place MultiObjectiveOptimizationFunction without cost_prototype" begin
+        function inplace_multi_obj_err!(cost, x, p)
+            cost[1] = sum(x .^ 2)
+            cost[2] = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return nothing
+        end
+        @test_throws ArgumentError MultiObjectiveOptimizationFunction(inplace_multi_obj_err!)
+    end
     sol = solve(prob, BBO_adaptive_de_rand_1_bin_radiuslimited(), callback = cb)
     # println(fitness_progress_history)
     @test !isempty(fitness_progress_history)

lib/OptimizationEvolutionary/test/runtests.jl

@@ -40,6 +40,44 @@ Random.seed!(1234)
         if state.iter % 10 == 0
             println(state.u)
         end
+
+    @testset "In-place Multi-Objective Optimization" begin
+        function inplace_multi_obj!(cost, x, p)
+            cost[1] = sum(x .^ 2)
+            cost[2] = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return nothing
+        end
+        u0 = [0.25, 0.25]
+        cost_prototype = zeros(2)
+        mof_inplace = MultiObjectiveOptimizationFunction(inplace_multi_obj!; cost_prototype=cost_prototype)
+        prob_inplace = OptimizationProblem(mof_inplace, u0)
+        sol_inplace = solve(prob_inplace, NSGA2())
+        @test sol_inplace ≠ nothing
+        @test length(sol_inplace.objective) == 2
+    end
+
+    @testset "Custom coalesce for Multi-Objective" begin
+        function multi_obj_tuple(x, p)
+            f1 = sum(x .^ 2)
+            f2 = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return (f1, f2)
+        end
+        coalesce_sum(cost, x, p) = sum(cost)
+        mof_coalesce = MultiObjectiveOptimizationFunction(multi_obj_tuple; coalesce=coalesce_sum)
+        prob_coalesce = OptimizationProblem(mof_coalesce, u0)
+        sol_coalesce = solve(prob_coalesce, NSGA2())
+        @test sol_coalesce ≠ nothing
+        @test mof_coalesce.coalesce([1.0, 2.0], [0.0, 0.0], nothing) == 3.0
+    end
+
+    @testset "Error if in-place MultiObjectiveOptimizationFunction without cost_prototype" begin
+        function inplace_multi_obj_err!(cost, x, p)
+            cost[1] = sum(x .^ 2)
+            cost[2] = sum(x .^ 2 .- 10 .* cos.(2π .* x) .+ 10)
+            return nothing
+        end
+        @test_throws ArgumentError MultiObjectiveOptimizationFunction(inplace_multi_obj_err!)
+    end
         return false
     end
     solve(prob, CMAES(μ = 40, λ = 100), callback = cb, maxiters = 100)