Manopt.jl

Manopt.jl is a package providing solvers for optimization problems defined on Riemannian manifolds. The implementation is based on ManifoldsBase.jl interface and can hence be used for all maniolds defined in Manifolds or any other manifold implemented using the interface.

Installation: OptimizationManopt.jl

To use the Optimization.jl interface to Manopt, install the OptimizationManopt package:

import Pkg;
Pkg.add("OptimizationManopt");

Methods

The following methods are available for the OptimizationManopt package:

• GradientDescentOptimizer: Corresponds to the gradient_descent method in Manopt.
• NelderMeadOptimizer : Corresponds to the NelderMead method in Manopt.
• ConjugateGradientDescentOptimizer: Corresponds to the conjugate_gradient_descent method in Manopt.
• ParticleSwarmOptimizer: Corresponds to the particle_swarm method in Manopt.
• QuasiNewtonOptimizer: Corresponds to the quasi_Newton method in Manopt.
• CMAESOptimizer: Corresponds to the cma_es method in Manopt.
• ConvexBundleOptimizer: Corresponds to the convex_bundle_method method in Manopt.
• FrankWolfeOptimizer: Corresponds to the FrankWolfe method in Manopt.

The common kwargs maxiters, maxtime and abstol are supported by all the optimizers. Solver specific kwargs from Manopt can be passed to the solve function or OptimizationProblem.

Examples

The Rosenbrock function on the Euclidean manifold can be optimized using the GradientDescentOptimizer as follows:

using Optimization, OptimizationManopt, Manifolds, LinearAlgebra
rosenbrock(x, p) = (p[1] - x[1])^2 + p[2] * (x[2] - x[1]^2)^2
x0 = zeros(2)
p = [1.0, 100.0]

R2 = Euclidean(2)

stepsize = Manopt.ArmijoLinesearch(R2)
opt = OptimizationManopt.GradientDescentOptimizer()

optf = OptimizationFunction(rosenbrock, Optimization.AutoZygote())

prob = OptimizationProblem(
    optf, x0, p; manifold = R2, stepsize = stepsize)

sol = Optimization.solve(prob, opt)

retcode: Failure
u: 2-element Vector{Float64}:
 0.8165527190561593
 0.6672818138967234

The box-constrained Karcher mean problem on the SPD manifold with the Frank-Wolfe algorithm can be solved as follows:

M = SymmetricPositiveDefinite(5)
m = 100
σ = 0.005
q = Matrix{Float64}(I, 5, 5) .+ 2.0
data2 = [exp(M, q, σ * rand(M; vector_at = q)) for i in 1:m]

f(x, p = nothing) = sum(distance(M, x, data2[i])^2 for i in 1:m)
optf = OptimizationFunction(f, Optimization.AutoZygote())
prob = OptimizationProblem(optf, data2[1]; manifold = M, maxiters = 1000)

function closed_form_solution!(M::SymmetricPositiveDefinite, q, L, U, p, X)
    # extract p^1/2 and p^{-1/2}
    (p_sqrt_inv, p_sqrt) = Manifolds.spd_sqrt_and_sqrt_inv(p)
    # Compute D & Q
    e2 = eigen(p_sqrt_inv * X * p_sqrt_inv) # decompose Sk  = QDQ'
    D = Diagonal(1.0 .* (e2.values .< 0))
    Q = e2.vectors

    Uprime = Q' * p_sqrt_inv * U * p_sqrt_inv * Q
    Lprime = Q' * p_sqrt_inv * L * p_sqrt_inv * Q
    P = cholesky(Hermitian(Uprime - Lprime))
    z = P.U' * D * P.U + Lprime
    copyto!(M, q, p_sqrt * Q * z * Q' * p_sqrt)
    return q
end
N = m
U = mean(data2)
L = inv(sum(1 / N * inv(matrix) for matrix in data2))

optf = OptimizationFunction(f, Optimization.AutoZygote())
prob = OptimizationProblem(optf, U; manifold = M, maxiters = 1000)

sol = Optimization.solve(
    prob, opt, sub_problem = (M, q, p, X) -> closed_form_solution!(M, q, L, U, p, X))

retcode: Failure
u: 5×5 Matrix{Float64}:
 2.99961  1.99957  1.99958  1.99922  1.99956
 1.99957  2.99953  1.99952  1.99913  1.9995
 1.99958  1.99952  2.99954  1.99914  1.99953
 1.99922  1.99913  1.99914  2.99867  1.99908
 1.99956  1.9995   1.99953  1.99908  2.99948

This example is based on the example in the Manopt and Weber and Sra'22.

The following example is adapted from the Rayleigh Quotient example in ManoptExamples.jl. We solve the Rayleigh quotient problem on the Sphere manifold:

using Optimization, OptimizationManopt
using Manifolds, LinearAlgebra
using Manopt

n = 1000
A = Symmetric(randn(n, n) / n)
manifold = Sphere(n - 1)

cost(x, p = nothing) = -x' * A * x
egrad(G, x, p = nothing) = (G .= -2 * A * x)

optf = OptimizationFunction(cost, grad = egrad)
x0 = rand(manifold)
prob = OptimizationProblem(optf, x0, manifold = manifold)

sol = solve(prob, GradientDescentOptimizer())

retcode: Failure
u: 1000-element Vector{Float64}:
  0.0031962538197439715
  0.01734655331376537
 -0.016887809626472013
  0.02960072973351771
 -0.006619316000918681
 -0.011289580433377779
 -0.05304819613716375
  0.044167718838723946
 -0.006421429589835838
 -0.03404482490840664
  ⋮
  0.03403219540802137
 -0.06293336520270104
  0.0746103077250663
 -0.0090803372680921
  0.02453855253911655
 -0.003119004570483792
 -0.018130919106824732
  0.019067454951662902
  0.019176241985761674

Let's check that this indeed corresponds to the minimum eigenvalue of the matrix A.

@show eigmin(A)
@show sol.objective

-0.06339392859867532