Adding support for new manifold to MANOPT

[Arrigo Benedetti]

Is it easy to add support for new manifold to MANOPT ? Is there any document/tutorial/example about how to do this? I would like to add support for the unitriangular matrix group.

Thanks

-Arrigo

[Nicolas Boumal]

Hello Arrigo,

Do you mean the set of square matrices with ones on the diagonal, zeros below the diagonal, and arbitrary real entries above the diagonal?

If so, since this is an affine space, it should be very easy to add support for this in manopt.

I suggest to take a look at this example:

https://github.com/NicolasBoumal/manopt/blob/master/manopt/manifolds/euclidean/symmetricfactory.m 

This is a "factory" (a file that adds support for a manifold) for the set of symmetric matrices of size n (the file also gives support for the manifold of tuples with k symmetric matrices each of size n).

Don't hesitate to let us know if anything is unclear, we'll be happy to help.

Best,

Nicolas

[Nicolas Boumal]

Ah, also: the various functions that you add to the structure M (such as M.inner, M.retr, M.rand etc) are described here:

https://www.manopt.org/tutorial.html#manifolds

(See after the list of available manifolds.)

[Arrigo Benedetti]

Nicolas,

thanks for your prompt reply. Yes, by unitriangular matrices I mean exactly what you wrote. I will look into the code for the symmetric matrix manifold as you suggested.

-Arrigo

[Nicolas Boumal]

If you get a good implementation for this factory that could be useful to others, feel free to contribute it to the toolbox; at any rate, we can help out here if you run into trouble.

--
http://www.manopt.org
---
You received this message because you are subscribed to the Google Groups "Manopt" group.
To unsubscribe from this group and stop receiving emails from it, send an email to manopttoolbo...@googlegroups.com.
Visit this group at https://groups.google.com/group/manopttoolbox.
To view this discussion on the web visit https://groups.google.com/d/msgid/manopttoolbox/b75f3b8d-720a-40ca-a162-48f54946362d%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Arrigo Benedetti]

Nicolas,

I just realized that I can reformulate my problem using an element of SO(n). My problem essentially consists in finding the maximal volume hyper-ellipsoid contained in a given convex polytope with constraints on the ratios of the ellipsoid semi-axes. Unfortunately this semi-axes constraint makes the problem nonconvex:

max log t

subject to

t *(a_i^T U diag(1,l_2,...,l_n) U^T a_i)^(1/2) + a_i^T c <= b_i,    i=1:m

where

U \in SO(n)

l_2,...,l_n are the known semiaxes ratios

a_i and b_i (i=1:m) are resp. the known normal vectors and distances to the origin of the hyperplanes defining the convex polytope

The ellipsoid is defined by

E = {x | x = c + t U diag(1,l_2,...,l_n) U^T w , norm(w) <= 1}

So now my question is: is manopt able to handle these kind of problems?

Thanks much

-Arrigo

[Nicolas Boumal]

Hello Arrigo,

From what I understood, you have a fixed ellipsoid, which you are allowed to (1) scale, and (2) rotate. The constraint is that the resulting ellipsoid must be included inside a given convex polytope, and the goal is to maximize the volume of the resulting ellipsoid.

Assuming the center of the ellipsoid is strictly included in the polytope, based on your latest message, we should be able to define a function like so:

For any rotation R in SO(d),

f(R) = the maximum scaling factor by which we can inflate the reference ellipsoid and still be included in the polytope.

Since the corresponding volume is an increasing function of f(R), it is equivalent to maximize f(R) or the volume associated to orientation R and scaling f(R).

The function f(R), I expect, will be nonsmooth, but might not be too hard to compute. I expect it should be possible to compute it exactly, but worst-case scenario, a bisection method would do the trick.

Then, to maximize f(R) over SO(d) in Manopt, you could use one of the derivative-free solvers such as pso or neldermead. These tend not to be particularly efficient, but since the cost function here is nonsmooth and getting its derivatives will be challenging, that's your best first bet. For rotation groups in low dimension (small d), this should be fine, depending on how fast things have to be.

By the way, this approach is fairly similar to this paper about bounding boxes:

Fast oriented bounding box optimization on the rotation group SO(3, R)

CHIA-TCHE CHANG, BASTIEN GORISSEN and SAMUEL MELCHIOR

http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.717.9566&rep=rep1&type=pdf 

Technical note: manopt minimizes always, so to maximize f(R), minimize -f(R).

(As a side note: I would expect that convexity already broke from the fact that you are maximizing (rather than minimizing) a volume.)

[Arrigo Benedetti]

Nicolas,

thanks for you comments! I think that in order to compute f(R) one just needs to solve a Linear Programming problem. By the way, I had not realized that the cost function is not smooth! As far as your comment about what breaks the convexity, I think it is imposing the ratios for the semi-axes or equivalently the eigenvalues of the Symmetric Positive Semidefinite matrix B that defines the ellipsoid. As shown in Boyd & Vandenberghe (Sec. 8.4.2) in the general case (no constraints) one can equivalently minimize log det B^-1 or maximize log det B.

I will try to use the Nelder & Mead solver in manopt and cvx to to compute f(R).

On a different note: are there any plans to port manopt to Julia?

Thanks

-Arrigo

[Nicolas Boumal]

Hello Arrigo,

 - I don't think you need linear programming to compute f(R). Computing f(R) involves finding one parameter: the optimal scaling factor. For scaling factors in [0, optimal scaling] the ellipsoid is inside the polytope. For scaling factors in (optimal scaling, +infty), it is not. And it is easy to test for each possible value of the scaling factor whether it is in the first range or the second range. So, you can start with an initial guess. If it's "acceptable" (the scaled ellipsoid is in the elliptope), then double your guess. Do so until you get a violation. In the end, this gives you an interval [a, b] such that a is too small, and b is too large. Now do bisection on this interval. I'm not saying it's the most efficient way, but it should be easy to implement.

 - Yes! Ronny Bergmann from TU Chemnitz is hard at work leading the development of a version of Manopt in Julia. I'm not sure what the timeline is for a public release. Feel free to contact him directly if you would like to test the current version.