Manopt.jl integration

docs/Project.toml

@@ -9,6 +9,8 @@ Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
 Ipopt_jll = "9cc047cb-c261-5740-88fc-0cf96f7bdcc7"
 IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e"
 Juniper = "2ddba703-00a4-53a7-87a5-e8b9971dde84"
+Manifolds = "1cead3c2-87b3-11e9-0ccd-23c62b72b94e"
+Manopt = "0fc0a36d-df90-57f3-8f93-d78a9fc72bb5"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 NLopt = "76087f3c-5699-56af-9a33-bf431cd00edd"
 Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"
@@ -17,6 +19,7 @@ OptimizationCMAEvolutionStrategy = "bd407f91-200f-4536-9381-e4ba712f53f8"
 OptimizationEvolutionary = "cb963754-43f6-435e-8d4b-99009ff27753"
 OptimizationGCMAES = "6f0a0517-dbc2-4a7a-8a20-99ae7f27e911"
 OptimizationMOI = "fd9f6733-72f4-499f-8506-86b2bdd0dea1"
+OptimizationManopt = "e57b7fff-7ee7-4550-b4f0-90e9476e9fb6"
 OptimizationMetaheuristics = "3aafef2f-86ae-4776-b337-85a36adf0b55"
 OptimizationMultistartOptimization = "e4316d97-8bbb-4fd3-a7d8-3851d2a72823"
 OptimizationNLopt = "4e6fcdb7-1186-4e1f-a706-475e75c168bb"

docs/pages.jl

@@ -23,6 +23,7 @@ pages = ["index.md",
         "Evolutionary.jl" => "optimization_packages/evolutionary.md",
         "Flux.jl" => "optimization_packages/flux.md",
         "GCMAES.jl" => "optimization_packages/gcmaes.md",
+        "Manopt.jl" => "optimization_packages/manopt.md",
         "MathOptInterface.jl" => "optimization_packages/mathoptinterface.md",
         "Metaheuristics.jl" => "optimization_packages/metaheuristics.md",
         "MultistartOptimization.jl" => "optimization_packages/multistartoptimization.md",

docs/src/examples/rosenbrock.md

@@ -158,6 +158,21 @@ using OptimizationBBO
 prob = Optimization.OptimizationProblem(rosenbrock, [0.0, 0.3], _p, lb = [-1.0, 0.2],
     ub = [0.8, 0.43])
 sol = solve(prob, BBO_adaptive_de_rand_1_bin()) # -1.0 ≤ x[1] ≤ 0.8, 0.2 ≤ x[2] ≤ 0.43
+
+## Riemannian optimization
+
+using Manopt, Manifolds, OptimizationManopt
+
+function rosenbrock_grad!(storage, x, p)
+    ForwardDiff.gradient!(storage, y -> rosenbrock(y, p), x)
+    project!(Manifolds.Sphere(1), storage, x, storage)
+end
+
+optprob = OptimizationFunction(rosenbrock; grad=rosenbrock_grad!)
+opt = OptimizationManopt.QuasiNewtonOptimizer(Manifolds.Sphere(1))
+x0 = [1.0, 0.0]
+prob = OptimizationProblem(optprob, x0, _p)
+sol = Optimization.solve(prob, opt)
 ```
 
 And this is only a small subset of what Optimization.jl has to offer!

docs/src/index.md

@@ -57,6 +57,7 @@ packages.
 | Evolutionary            | ❌                    | ❌                   | ❌                     | ✅               | ❌                 | ✅                    | 🟡                    |
 | Flux                    | ✅                    | ❌                   | ❌                     | ❌               | ❌                 | ❌                    | ❌                    |
 | GCMAES                  | ❌                    | ❌                   | ❌                     | ✅               | ❌                 | ✅                    | ❌                    |
+| Manopt                  | ✅                    |                   🟡 | ✅                     | ✅               | ✅                 | ✅                    | ❌                    |
 | MathOptInterface        | ✅                    | ✅                   | ✅                     | ✅               | ✅                 | ✅                    | 🟡                    |
 | MultistartOptimization  | ❌                    | ❌                   | ❌                     | ✅               | ❌                 | ✅                    | ❌                    |
 | Metaheuristics          | ❌                    | ❌                   | ❌                     | ✅               | ❌                 | ✅                    | 🟡                    |

docs/src/optimization_packages/manopt.md

@@ -0,0 +1,81 @@
+# [Manopt.jl](@id manopt)
+[`Manopt`](https://github.com/JuliaManifolds/Manopt.jl) is Julia package implementing various algorithm focusing on manifold-based constraints. When all or some of the constraints are in the form of Riemannian manifolds, optimization can be often performed more efficiently than with generic methods handling equality or inequality constraints.
+
+See [`Manifolds`](https://github.com/JuliaManifolds/Manifolds.jl) for a library of manifolds that can be used as constraints.
+
+## Installation: OptimizationManopt.jl
+
+To use this package, install the `OptimizationManopt` package:
+
+```julia
+import Pkg; Pkg.add("OptimizationManopt")
+```
+
+## Methods
+
+`Manopt.jl` algorithms can be accessed via Optimization.jl using one of the following optimizers:
+
+- `ConjugateGradientDescentOptimizer`
+- `GradientDescentOptimizer`
+- `NelderMeadOptimizer`
+- `ParticleSwarmOptimizer`
+- `QuasiNewtonOptimizer`
+
+For a more extensive documentation of all the algorithms and options please consult the [`Documentation`](https://manoptjl.org/stable/).
+
+## Local Optimizer
+
+### Derivative-Free
+
+The Nelder-Mead optimizer can be used for local derivative-free optimization on manifolds.
+
+```@example Manopt1
+using Optimization, OptimizationManopt, Manifolds
+rosenbrock(x, p) =  (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2
+x0 = [0.0, 1.0]
+p  = [1.0,100.0]
+f = OptimizationFunction(rosenbrock)
+opt = OptimizationManopt.NelderMeadOptimizer(Sphere(1))
+prob = Optimization.OptimizationProblem(f, x0, p)
+sol = solve(prob, opt)
+```
+
+### Gradient-Based
+
+Manopt offers gradient descent, conjugate gradient descent and quasi-Newton solvers for local gradient-based optimization.
+
+Note that the returned gradient needs to be Riemannian, see for example [https://manoptjl.org/stable/functions/gradients/](https://manoptjl.org/stable/functions/gradients/).
+Note that one way to obtain a Riemannian gradient is by [projection and (optional) Riesz representer change](https://juliamanifolds.github.io/Manifolds.jl/latest/features/differentiation.html#Manifolds.RiemannianProjectionBackend).
+
+```@example Manopt2
+using Optimization, OptimizationManopt, Manifolds
+rosenbrock(x, p) =  (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2
+function rosenbrock_grad!(storage, x, p)
+    storage[1] = -2.0 * (p[1] - x[1]) - 4.0 * p[2] * (x[2] - x[1]^2) * x[1]
+    storage[2] = 2.0 * p[2] * (x[2] - x[1]^2)
+    project!(Sphere(1), storage, x, storage)
+end
+x0 = [0.0, 1.0]
+p  = [1.0,100.0]
+f = OptimizationFunction(rosenbrock; grad = rosenbrock_grad!)
+opt = OptimizationManopt.GradientDescentOptimizer(Sphere(1))
+prob = Optimization.OptimizationProblem(f, x0, p)
+sol = solve(prob, opt)
+```
+
+## Global Optimizer
+
+### Without Constraint Equations
+
+The particle swarm optimizer can be used for global optimization on a manifold without constraint equations. It can be especially useful on compact manifolds such as spheres or orthogonal matrices.
+
+```@example Manopt3
+using Optimization, OptimizationManopt, Manifolds
+rosenbrock(x, p) =  (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2
+x0 = [0.0, 1.0]
+p  = [1.0,100.0]
+f = OptimizationFunction(rosenbrock)
+opt = OptimizationManopt.ParticleSwarmOptimizer(Sphere(1))
+prob = Optimization.OptimizationProblem(f, x0, p)
+sol = solve(prob, opt)
+```

docs/src/tutorials/constraints.md

@@ -4,7 +4,7 @@ Multiple optimization packages available with the MathOptInterface and Optim's `
 Optimization.jl provides a simple interface to define the constraint as a Julia function and then specify the bounds for the output
 in `OptimizationFunction` to indicate if it's an equality or inequality constraint.
 
-Let's define the rosenbrock function as our objective function and consider the below inequalities as our constraints.
+Let's define the Rosenbrock function as our objective function and consider the below inequalities as our constraints.
 
 ```math
 \begin{aligned}
@@ -99,3 +99,42 @@ res = zeros(2)
 cons(res, sol.u, _p)
 println(res)
 ```
+
+## Constraints as Riemannian manifolds
+
+Certain constraints can be efficiently handled using Riemannian optimization methods. This is the case when moving the solution on a manifold can be done quickly. See [here](https://juliamanifolds.github.io/Manifolds.jl/latest/index.html) for a (non-exhaustive) list of such manifolds, most prominent of them being spheres, hyperbolic spaces, Stiefel and Grassmann manifolds, symmetric positive definite matrices and various Lie groups.
+
+Let's for example solve the Rosenbrock function optimization problem with just the spherical constraint. Note that the constraint isn't passed to `OptimizationFunction` but to the optimization method instead. Here we will use a [quasi-Newton optimizer based on the BFGS algorithm](https://manoptjl.org/stable/solvers/quasi_Newton/).
+
+```@example manopt
+using Optimization, Manopt, OptimizationManopt, Manifolds
+
+rosenbrock(x, p) = (p[1] - x[1])^2 + p[2] * (x[2] - x[1]^2)^2
+x0 = [0.0, 1.0]
+_p = [1.0, 100.0]
+
+function rosenbrock_grad!(storage, x, p)
+    # the first part can be computed using AD tools
+    storage[1] = -2.0 * (p[1] - x[1]) - 4.0 * p[2] * (x[2] - x[1]^2) * x[1]
+    storage[2] = 2.0 * p[2] * (x[2] - x[1]^2)
+    # projection is needed because Riemannian optimizers expect
+    # Riemannian gradients instead of Euclidean ones.
+    project!(Sphere(1), storage, x, storage)
+end
+
+optprob = OptimizationFunction(rosenbrock; grad=rosenbrock_grad!)
+opt = OptimizationManopt.QuasiNewtonOptimizer(Sphere(1))
+prob = OptimizationProblem(optprob, x0, _p)
+sol = Optimization.solve(prob, opt)
+```
+
+Note that currently `AutoForwardDiff` can't correctly compute the required Riemannian gradient for optimization. Riemannian optimizers require Riemannian gradients while `ForwardDiff.jl` returns normal Euclidean ones. Conversion from Euclidean to Riemannian gradients can be performed using the [`project`](https://juliamanifolds.github.io/ManifoldsBase.jl/stable/projections.html#Projections) function and (for certain manifolds) [`change_representer`](https://juliamanifolds.github.io/Manifolds.jl/stable/manifolds/metric.html#Manifolds.change_representer-Tuple{AbstractManifold,%20AbstractMetric,%20Any,%20Any}).
+
+Note that the constraint is correctly preserved and the convergence is quite fast.
+
+```@example manopt
+println(norm(sol.u))
+println(sol.original.stop.reason)
+```
+
+It is possible to mix Riemannian and equation-based constraints but it is currently a topic of active research. Manopt.jl offers solvers for such problems but they are not accessible via the Optimization.jl interface yet.

lib/OptimizationManopt/LICENSE

@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) 2022 Julia Manifolds
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.

lib/OptimizationManopt/Project.toml

@@ -0,0 +1,10 @@
+name = "OptimizationManopt"
+uuid = "e57b7fff-7ee7-4550-b4f0-90e9476e9fb6"
+authors = ["Mateusz Baran <[email protected]>"]
+version = "0.1.0"
+
+[deps]
+Manifolds = "1cead3c2-87b3-11e9-0ccd-23c62b72b94e"
+ManifoldsBase = "3362f125-f0bb-47a3-aa74-596ffd7ef2fb"
+Manopt = "0fc0a36d-df90-57f3-8f93-d78a9fc72bb5"
+Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"