Manopt.jl

Manopt.jl is a package with implementations of a variety of optimization solvers on manifolds supported by Manifolds.

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}:
 3.00129  2.00134  2.00149  2.00146  2.00149
 2.00134  3.00151  2.0016   2.00156  2.00147
 2.00149  2.0016   3.00165  2.00164  2.00163
 2.00146  2.00156  2.00164  3.00157  2.00163
 2.00149  2.00147  2.00163  2.00163  3.00155

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.014896160311621047
 -0.001497482509883409
 -0.00313087508504342
  0.026490267884331647
 -0.0008612539385312455
 -0.046585705164980874
 -0.049234311732548004
  0.023933947754123983
  0.019545978556177773
 -0.06311516665107293
  ⋮
 -0.028454214726346662
  0.027732939452908736
  0.044220828165509735
 -0.026074843309303125
  0.04950109816608916
  0.002267480827377856
 -0.00351459583146516
  0.016076884193200816
  0.0038319767558302794

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

@show eigmin(A)
@show sol.objective

-0.06306679120698035