fixes tests

SciML · Oct 21, 2022 · 9a8ef78 · 9a8ef78
1 parent 90a033a
commit 9a8ef78
Show file tree

Hide file tree

Showing 4 changed files with 137 additions and 100 deletions.
diff --git a/docs/Project.toml b/docs/Project.toml
@@ -3,11 +3,12 @@ BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
 Optim = "429524aa-4258-5aef-a3af-852621145aeb"
+Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"
+OptimizationOptimJL = "36348300-93cb-4f02-beb5-3c3902f8871e"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"

diff --git a/docs/src/tutorials/optimization.md b/docs/src/tutorials/optimization.md
@@ -42,7 +42,7 @@ p = [
 ]
 
 prob = OptimizationProblem(sys, u0, p, grad=true, hess=true)
-solve(prob, Newton())
+solve(prob, GradientDescent())
 ```
 
 ## Rosenbrock Function with Constraints
@@ -63,8 +63,8 @@ u0 = [
     x => 1.0
     y => 2.0
 ]
-prob = OptimizationProblem(sys, u0, p, grad=true, hess=true)
-solve(prob, IpNewton())
+prob = OptimizationProblem(sys, u0, grad=true, hess=true)
+solve(prob, IPNewton())
 ```
 
 A visualization of the objective function and the inequality constraint is depicted below.
@@ -76,10 +76,10 @@ contour(x, y, (x,y) -> (1 - x)^2 + 100 * (y - x^2)^2, fill=true, color=:viridis,
 contour!(x, y, (x, y) -> x^2 + y^2, levels=[1], color=:lightblue, line=4)
 ```
 
-### Explenation
+### Explanation
 Equality and inequality constraints can be added to the `OptimizationSystem`. An equality constraint can be specified via and `Equation`, e.g., `x^2 + y^2 ~ 1`. While inequality constraints via an `Inequality`, e.g., `x^2 + y^2 ≲ 1`. The syntax is here `\lesssim` and `\gtrsim`.
 
-## Nexted Systems
+## Nested Systems
 Needs more text but it's super cool and auto-parallelizes and sparsifies too.
 Plus you can hierarchically nest systems to have it generate huge
 optimization problems.
diff --git a/src/systems/optimization/optimizationsystem.jl b/src/systems/optimization/optimizationsystem.jl
@@ -223,14 +223,25 @@ function DiffEqBase.OptimizationProblem{iip}(sys::OptimizationSystem, u0map,
     ps = parameters(sys)
     cstr = constraints(sys)
 
+    lb = first.(getbounds.(dvs))
+    ub = last.(getbounds.(dvs))
+    int = isintegervar.(dvs) .| isbinaryvar.(dvs)
+    lb[isbinaryvar.(dvs)] .= 0
+    ub[isbinaryvar.(dvs)] .= 1
+
     defs = defaults(sys)
     defs = mergedefaults(defs, parammap, ps)
     defs = mergedefaults(defs, u0map, dvs)
 
     u0 = varmap_to_vars(u0map, dvs; defaults = defs, tofloat = false)
     p = varmap_to_vars(parammap, ps; defaults = defs, tofloat = false, use_union)
-    lb = varmap_to_vars(lb, dvs; check = false, tofloat = false, use_union)
-    ub = varmap_to_vars(ub, dvs; check = false, tofloat = false, use_union)
+    lb = varmap_to_vars(dvs .=> lb, dvs; defaults = defs, tofloat = false, use_union)
+    ub = varmap_to_vars(dvs .=> ub, dvs; defaults = defs, tofloat = false, use_union)
+
+    if all(lb .== -Inf) && all(ub .== Inf)
+        lb = nothing
+        ub = nothing
+    end
 
     f = generate_function(sys, checkbounds = checkbounds, linenumbers = linenumbers,
                           expression = Val{false})
@@ -270,10 +281,6 @@ function DiffEqBase.OptimizationProblem{iip}(sys::OptimizationSystem, u0map,
         hess_prototype = nothing
     end
 
-    # TODO: handle integrality constraints
-    lb = first.(getbounds.(dvs)) # TODO: handle case where bounds are parameters
-    ub = last.(getbounds.(dvs))
-
     if length(cstr) > 0
         @named cons_sys = ConstraintsSystem(cstr, dvs, ps)
         cons, lcons, ucons = generate_function(cons_sys, checkbounds = checkbounds,
@@ -308,7 +315,7 @@ function DiffEqBase.OptimizationProblem{iip}(sys::OptimizationSystem, u0map,
                                                   cons_hess_prototype = cons_hess_prototype,
                                                   expr = obj_expr,
                                                   cons_expr = cons_expr)
-        OptimizationProblem{iip}(_f, u0, p; lb = lb, ub = ub,
+        OptimizationProblem{iip}(_f, u0, p; lb = lb, ub = ub, int = int,
                                  lcons = lcons, ucons = ucons, kwargs...)
     else
         _f = DiffEqBase.OptimizationFunction{iip}(f,
@@ -320,7 +327,7 @@ function DiffEqBase.OptimizationProblem{iip}(sys::OptimizationSystem, u0map,
                                                   paramsyms = Symbol.(parameters(sys)),
                                                   hess_prototype = hess_prototype,
                                                   expr = obj_expr)
-        OptimizationProblem{iip}(_f, u0, p; lb = lb, ub = ub,
+        OptimizationProblem{iip}(_f, u0, p; lb = lb, ub = ub, int = int,
                                  kwargs...)
     end
 end
@@ -394,14 +401,25 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
         hess_prototype = nothing
     end
 
+    lb = first.(getbounds.(dvs))
+    ub = last.(getbounds.(dvs))
+    int = isintegervar.(dvs) .| isbinaryvar.(dvs)
+    lb[isbinaryvar.(dvs)] .= 0
+    ub[isbinaryvar.(dvs)] .= 1
+
     defs = defaults(sys)
     defs = mergedefaults(defs, parammap, ps)
     defs = mergedefaults(defs, u0map, dvs)
 
     u0 = varmap_to_vars(u0map, dvs; defaults = defs, tofloat = false)
     p = varmap_to_vars(parammap, ps; defaults = defs, tofloat = false, use_union)
-    lb = varmap_to_vars(lb, dvs; check = false, tofloat = false, use_union)
-    ub = varmap_to_vars(ub, dvs; check = false, tofloat = false, use_union)
+    lb = varmap_to_vars(dvs .=> lb, dvs; defaults = defs, tofloat = false, use_union)
+    ub = varmap_to_vars(dvs .=> ub, dvs; defaults = defs, tofloat = false, use_union)
+
+    if all(lb .== -Inf) && all(ub .== Inf)
+        lb = nothing
+        ub = nothing
+    end
 
     obj_expr = toexpr(objective(sys))
     pairs_arr = if p isa SciMLBase.NullParameters
@@ -412,10 +430,6 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
     end
     rep_pars_vals!(obj_expr, pairs_arr)
 
-    # TODO: handle integrality constraints
-    lb = first.(getbounds.(dvs)) # TODO: handle case where bounds are parameters
-    ub = last.(getbounds.(dvs))
-
     if length(cstr) > 0
         @named cons_sys = ConstraintsSystem(cstr, dvs, ps)
         cons, lcons, ucons = generate_function(cons_sys, checkbounds = checkbounds,
@@ -443,6 +457,7 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
             hess = $_hess
             lb = $lb
             ub = $ub
+            int = $int
             cons = $cons[1]
             lbcons = $lbcons
             ubcons = $ubcons
@@ -463,7 +478,7 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
                                            cons_hess_prototype = cons_hess_prototype,
                                            expr = obj_expr,
                                            cons_expr = cons_expr)
-            OptimizationProblem{$iip}(_f, u0, p; lb = lb, ub = ub, lcons = lcons,
+            OptimizationProblem{$iip}(_f, u0, p; lb = lb, ub = ub, int = int, lcons = lcons,
                                       ucons = ucons, kwargs...)
         end
     else
@@ -475,6 +490,7 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
             hess = $_hess
             lb = $lb
             ub = $ub
+            int = $int
             syms = $(Symbol.(states(sys)))
             paramsyms = $(Symbol.(parameters(sys)))
             _f = OptimizationFunction{iip}(f, SciMLBase.NoAD();
@@ -484,7 +500,7 @@ function OptimizationProblemExpr{iip}(sys::OptimizationSystem, u0,
                                            paramsyms = paramsyms,
                                            hess_prototype = hess_prototype,
                                            expr = obj_expr)
-            OptimizationProblem{$iip}(_f, u0, p; lb = lb, ub = ub, kwargs...)
+            OptimizationProblem{$iip}(_f, u0, p; lb = lb, ub = ub, int = int, kwargs...)
         end
     end
 end