[ANN] Manopt.jl

The new version, Manopt.jl 0.2.1 is now using ManifoldsBase.jl as its backend, so all manifolds from Manifolds.jl, see also announcement here. It also has a new home at JuliaManifolds, our new org on GitHub. This for example also resolves the discussions on parallel and other vector transports, since those are now covered in the general manifold case.

11 Likes

The newest version 0.3.0 of Manopt.jl is a major rework, since one can now also provide mutating functions gradF!(M, X, p) for the gradient, avoiding allocations of tangent vectors every iteration. A new tutorial https://manoptjl.org/stable/tutorials/Benchmark.html illustrates the new feature, that is available for all functions a problem can contain despite the cost itself, i.e. for gradients, proximal maps, hessians, differentials, adjoint differentials. All functions form the library were also updated accordingly.

8 Likes

If you are working with Manopt.jl you can now cite the new JOSS paper at Journal of Open Source Software: Manopt.jl: Optimization on Manifolds in Julia

10 Likes

The newest release of Manopt.jl 0.3.47 covers constrained optimisation tasks on Manifolds as well, see

as well as the tutorial at Do constrained Optimization · Manopt.jl.

6 Likes

It’s been more than 1.5 years since the last breaking change in Manopt.jl!
We introduced several algorithms in that time, but now we decided it was time for some major changes. We introduce

Manopt.jl 0.4

The main new features are

  • A type for the objective that collects the cost and all ist dependent functions like the gradient, proximal maps or the Hessian
  • A first simple cache for the objective, see Cache Objectives. We plan to extend these caching capabilities soon by a general cache.
  • random points and tangent vectors are now based on the definition in ManifoldsBase.jl making the use of random points even simpler and more variants from Manifolds.jl available
  • the new objective allows to introduce the CostGrad structure, when you have one function (e.g. from AD tools) that simultaneously evaluates the cost and the gradient
  • instead of definiing own differential basic functions, we now depend on ManifoldDiff.jl, which provides more flexibility in using differentiation tools now and in the future.
  • We unified the naming scheme with Manifolds.jl

Recently we also introduced

If you want to get started, check out the Optimize! tutorial, also available as Pluto notebook.

9 Likes

The most recent version 0.4.8 introduces a comprehensive display of the resulting solver state. To illustrate this consider a small example on computing the Riemannian center of mass

using Manopt, Manifolds
M = Sphere(2)
n=100;
p = 1 / sqrt(2) * [1.0, 0.0, 1.0];
data = [exp(M, p,  π / 8 * rand(M; vector_at=p)) for _ in 1:n];
f(M, p) = sum(1 / (2 * n) * distance.(Ref(M), Ref(p), data) .^ 2);
grad_f(M, p) = sum(1 / n * grad_distance.(Ref(M), data, Ref(p)));
p_opt = gradient_descent(M, f, grad_f, data[1])

The solver call, that is the p_opt you obtain is (just) the minimiser the gradient descent obtained. So if you just need that, this is as before.
Also before you could ask for the whole solver state at the end of the solver run using

state = gradient_descent(M, f, grad_f, data[1]; return_state=true)

But now these state structures together with a lot of internal structures have new/updated/nicer show methods and an even shorter status_summary() function such that this state gets displayed on REPL as

# Solver state for `Manopt.jl`s Gradient Descent
After 75 iterations

## Parameters
* retraction method: ExponentialRetraction()

## Stepsize
ArmijoLineseach() with keyword parameters
  * initial_stepsize = 1.0
  * retraction_method = ExponentialRetraction()
  * contraction_factor = 0.95
  * sufficient_decrease = 0.1
  * linesearch_stopsize = 0.0

## Stopping Criterion
Stop When _one_ of the following are fulfilled:
    Max Iteration 200:  not reached
    |Δf| < 1.0e-9: reached
Overall: reached
This indicates convergence: Yes

also including a heuristic whether the stopping criterion that did “fire” (indicate to stop) is one that indicates convergence. Note that you can easily design your own stopping criteria, since every solver in Manopt.jl has the stopping_criterion= keywords and stopping criteria can be combined using | and & for example if you wand both a small change in the iterate and a small gradient – or a number of iterations as a safeguard you can do

stopping_criterion = StopAfterIteration(200) | (
     StopWhenGradientNormLess(1e-9) & StopWhenChangeLess(1e-9)
)

Finally, the examples of this package will get their own “home” and have already “moved out” – they will be revised and improved in performance and then be added to ManoptExamples.jl, see this announcement.

8 Likes

We have two new algorithms in the newest version

Manopt.jl 0.4.12

namely

  • the Difference of Convex Algorithm (DCA)
  • the Difference of Convex Proximal Point Algorithm (DCPPA)

to solve problems that can be written as a difference of two convex functions (DC problems) – of course those phrased on Riemannian manifolds!
See Difference of Convex · Manopt.jl for details on both these algorithms

8 Likes

We are making small further steps in small features.

For example this week we implemented a new check in Manopt 0.4.18 so that you can numerically check you Implementation of the Hessian, see
https://manoptjl.org/stable/helpers/checks/

for Details (and of course to look at the theoretical source by Nicolas Boumal)

Overall, the new 0.4 from January this year seems to be the interface we probably stick with, it seems nice and easy to use. We of course still have some things to do Issues · JuliaManifolds/Manopt.jl · GitHub – feel free to add your ideas as well, but I feel we are on a good track to a 1.0.
The main feature I am missing is maybe even more something for ManifoldDiff, namely conversion of Euclidean to Riemannian Hessians – just if someone wants to join into that adventure :slight_smile:

2 Likes

I am using Manopt extensively, and find that the Steinhaug-Toint trust region subproblem solver is not as efficient as using the Lanczos method. I have large Hessian matrices, and even when I fully calculate them, the Lanczos method is much better than the supposedly sparse Steinhaug-Toint. Is there a plan to add different implementations, or the focus is only on the interface?

1 Like

Oh, nice to hear!

What do you mean by “only on the interface”? We currently provide more than 20 solvers, so I would say that that might be a sign, that we focus on both in equal manners.

One challenge here is, that all interface design and 90% of the solvers are implementation by just myself. The remaining 10% were students of mine (3-4 master theses) – currently for example an ARC solver that should be faster than TR.

for TR and Steinhaug-Toint, that was also a student project and some while ago – definetly not the students fault but also due to my experience in speed, the solver might not be that fast.
And sure, if you compare to ALM or the more recent DC solvers, the trust region one still misses the flexibility that currently the suit solver is not so easy exchangeable. But it would be great to have, for sure to

  1. make that more flexible that the subsolver can be exchanged (I added an issue here since this hopefully is more like a technical thing)
  2. implement different solves for the sub problem in trust region (I have no clue about these).

But I have to honestly say that for the next year or so, I do not see much time in my schedule (maybe I can do that with a student project).
Keep in mind that most of that is done in my free time – at least the algorithms I do not do research on myself.

Also – contributions are of course always welcome – I can definelty help with step 1 – but for different sub solvers I would also have to read on them and their efficiency.

edit: Or to be precise – this mechanism of specifying sub solvers was developed last autumn – while the trust region solver with its sub solver was written four years ago. So it can definetly be improved both in flexibility for sub solvers and in speed itself – it is just a time issue on my side. So feel free to help :slight_smile:

edit2: If you have references and can describe what you need – for best of cases with manifold references, feel free to open that as an issue

2 Likes

Thanks for the quick reply. Yes, time is a major problem, I do not have much either. Would it be possible to just plug in this trust region solver instead of the current one?

It implements a matrix-free method.

In my own package (see GitHub - rs1909/FMA: Foliations, Manifolds and Autoencoders) I need to use a huge productmanifold. I then use a factor-wise optimisation (a.k.a. batch coordinate descent), where I implemented a trust region solver myself using TRS.jl. But I calculate the full Hessian for each matrix manifold component. It would be nice if I did not have to do that.

Well, time-wise, I already implemented quite a few solvers on manifolds just out of curiosity, but since I do not use all of them – for example I often use L-BFGS but not really often TR – I also not often check and improve them.

Product manifolds we cover quite nicely in Product manifold · Manifolds.jl where we came up with a very nice and efficient representation/storage.

Concerning the linked solver – one challenge there is that the Manopt.jl TR-subsolver is a subsolver in a tangent space. While ManifoldsBase.jl does provide most methods from vector spaces for our tangent spaces – and TRS.jl seems to be a Euclidean subsolver? That would at least require some kind of wrapper – for security but also to match the interface.

edit: Oh, no its even worse, the radius constraint in the link you said, has of course to be the radius with respect to the norm in the tangent space (induced by the Riemannian inner product) – usually we never set up the full metric tensor that we would need to pass to your ellipsoidal norms?
And of course for a wrapper (and a weak dependency) the package would have to be registered.

Sure such a wrapper would be possible after step 1 mentioned above – one idea of the sub solvers Manopt.jl currently has in place is that for example the storage for the iterate in subproblems can be reused in the next subsolver call (no new allocations necessary) and all Manifolds/Manopt functions usually work in-place (also Hessian Computations). One thing you could consider is maybe approximate hessians (cf. Trust-Regions Solver · Manopt.jl)? But sure instead of those if the Hessian is Sparse, that might also be beneficial.

So, with a bit of caveat (Eclidean vs Tangent space) that should be doable – also the wrapper idea for other solvers is nice. It just remains a time issue – so step 2. Step 1 is not an opened issue and I try to work on these when I find time – the check of Hessians was one of these I had in mind to do when finding time :slight_smile:

edit 2: The story with the inner product x‘Cx is even more complicated since on general manifold it of course can also change depending on the point you are at; and since (again) we usually work chart-free, we even never explicitly write out C nor do we have x in a form of a vector of length equal to the dimension of the manifold. So a wrapper to your solver might indeed be complicated.

Maybe I’m missing something but isn’t the current trust-region subsolver in Manopt.jl essentially matrix-free? At least as long as you provide your own Hess_f that exploits sparsity, which seems like the way to go here.

Correct, the current Steinhaug-Toint does not build a whole matrix and if you implement Hessian efficiently this might be a way around.

I still read the remark as – we could also do Lanczos – which is probably an idea though I do not recall the Riemannian paper about that, but ARC does Lanczos and ARC is a more advanced TR.

So you are right with a (probably also inlace) efficient implementation of the Hessian already ST can be improved (without changing Manopt.jl) – still making the solver more flexible and providing a few might be nice.

Maybe @fastwave could open a thread about the implementation of the Hessian and we could discuss there (if you are interested, of course).

I am sorry that my comments were all over the place. The Steinhaug-Toint is definitely matrix free, but it is only approximate and also inefficient. It is good when an interior point is found but tends to be bad if the boundary is reached. A better method is in Adachi et al, Solving the trust regions subproblem… I am looking into how to implement it in Manopt…

Unfortunately I did not get the reference ARC. Could @kellertuer please clarify? If that could help with the implementation I would be very happy. :slight_smile:

I looked if someone already worked on a similar problem, and it turns out that there is an interesting algorithm: A limited-memory Riemannian symmetric rank-one trust-region method with an efficient algorithm for its subproblem - ScienceDirect . I haven’t checked all details but they appear to solve their subproblem in a basis. You could then use any Euclidean solver, probably? I can help hook it up to Manopt.jl :slightly_smiling_face: Manifolds.jl luckily has bases.

As I said for the method you reference, an interesting point is, whether this can be done by just calling inner (instead of the matrix C you have in your docs); but I will try to work on the generalisation (to more sub solvers) soon then :slight_smile:

Concerning ARC – sorry, I was imprecise there. It refers to the Adaptive Regulasirzation with cubics, which on manifolds is this paper [1806.00065] Adaptive regularization with cubics on manifolds – and as far as I know the default algorithm used in Manopt/Matlab. As I said, this is work in progress – and might also help us since the subproblem is a more complicated one than TR and we do it with Lanczos, so maybe one can adapt (or directly use) that Lanczos also for TR.

At least the symmetric rank one approximation of the Hessian is available (both as an update in quasi Newton as well as a hessian approximation cf. Trust-Regions Solver · Manopt.jl ); but sure other tr sub solvers remain to be made possible as mentioned.

I’ve made a draft PR: Add customizable subsolver to trust_regions by mateuszbaran · Pull Request #241 · JuliaManifolds/Manopt.jl · GitHub :slight_smile: . @fastwave , this kind of extensibility should be enough for you to replace the subsolver I think.

the new version was just registered (so is available later today). For a different subsolver, you need to mainly mimic what Steinhaug-Toiint does

  • define a new solver state
  • implement initialize_solver!(mp::AbstractManoptProblem, s::MyNewState) (if theree is something to do, there is a default “do-nothing” available)
  • implement step_solver!(mp::AbstractManoptProblem, s::MyNewState)
  • implement get_solver_result(s:: MyNewState) (if the results is not stored in s.p)

(And maybe from a moderator side the discussion starting from the “can we have more subsolvers” could be its own topic, since we are spamming a bit this announcement thread – as super interesting as this topic is :thumbsup:)