# MFOPT - Optimization over Manifolds

Download the most current release: Version 1.0, 15/10/2013

This package contains the *MFOPT* MATLAB library that was used for the computations in the original ICCV 2013 paper.

**Features:**
MFOPT can be used to approximately solve variational problems with total variation-base regularization and arbitrary data terms for functions that are constrained to values in a manifold. The library includes implementations for

- Euclidean spaces R^n,
- n-dimensional unit spheres S^n,
- the Moebius band,
- the three-dimensional rotation group SO(3).

Other manifolds can be easily added due to the modular structure.

**Requirements:** 64-bit version of MATLAB 2013a running under 64-bit Windows 7 or later.

**Limitations:**

- Figures 3,4,6 in the original paper were generated using a different code base and are not included.
- To clearly see the discretization effects as in Fig. 2, the problems have to be solved to a high accuracy. Low accuracy solutions that are not fully converged can actually look better (although they are in fact incorrect solutions). For the paper, a commercial solver was used to obtain high accuracy. The package contains a version that uses the simpler primal-dual method but gives much lower accuracy.

**License:**
You may freely use and modify the source code for non-commercial research, as long as a reference is included in related publications. Please cite:

```
J. Lellmann, E. Strekalovskiy, S. Koetter, and D. Cremers: Total Variation Regularization for Functions with Values in a Manifold. In: International Conference of Computer Vision (ICCV), 2013. In press.
```

## Quick Start

- Download and unzip to a folder of your choice
- run any of the
`paper_*`

scripts from within MATLAB to solve the optimization problems, or the`vis_*`

scripts to visualize the results and generate the figures. A good starting point is`paper_moebius.m`

. - Stopping criteria are set to high accuracy by default and will take a long time to compute. Computation time can be reduced by lowering the maximum iteration count
`term_maxiter`

or increasing the target accuracy`term_relgap`

.

## Overview

- The
`paper_*`

scripts read the input data from the`input`

folder, and solve the corresponding optimization problems. - All output is stored in separate folders within the
`results`

folder, named by the date and time the script was run. - The
`results_important`

folder contains copies of the exact results shown in the paper. - The
`vis_*`

scripts use the data in the`results_important`

folder to generate the figures. The figures are stored in the`results`

folder. - The actual minimization happens in
`solve_manifold()`

. If`use_gpu = true`

(default), the actual iterations are carried out by the MEX function`mf_iterate`

. Despite the field name, the open implementation is currently CPU only, so no GPU is required.

## Adding your own manifold

All information about the manifold is passed to the solver in form of a `struct`

that is created by one of the `manifold_*`

functions. In order to add your own manifold, the best way is to modify one of the `manifold_*`

functions.

The manifold structure has the following fields:

`m.mdims`

- struct holding the manifold dimensions and various constants:`m.mdims.dmanifold`

- dimension of the manifold: n for R^n, m for the m-dimensional sphere S^m, 2 for the Moebius band, and so on.`m.mdims.gradients`

- number of “evaluation points” at which the gradient is evaluated. Denoted by`m`

in the paper.`m.mdims.points`

- size of the neighbourhood N_j in (10). Each gradient computation in one of the evaluation points uses the value of u at`m.mdims.points`

points. This is used together with`m.P`

(see below).`m.mdims.vecs`

- the number of vectors (labels) used to discretize the manifold. Denoted by`l`

in the paper.

`m.v`

- the vectors/points on the manifold, embedded into some R^n, that correspond to the labels. Denoted by z^k in the paper.`m.b`

- the weights`b`

used to integrate functions over the manifold as in (9).`m.P`

- defines the neighbourhood N_j for the evaluation points as in (10). For the j-th evaluation point y^j, the neighbourhood consists of all`x^{i+1}`

with`i`

in`P(:,j)`

.`P`

must be zero-based and of`int64`

type.`mf_iterate`

contains a hard-coded limit which requires size(P,1) ⇐ 25 if`use_gpu = true`

is used.`m.A,m.B`

- These are the matrices A and B from equation (14) in the paper. They define the way gradients are computed on the manifold.`m.dist(v1,v2)`

- function handle, computes the geodesic distance between a pair of points x,y on the manifold. This uses the same vector representation/embedding into R^n as`m.v`

.`m.mean(z,u)`

- function handle, computes the Fréchet mean (22) for the vectors`z`

and weights`u`

.`m.plot()`

- function handle, sets up a (2D or 3D) figure and axis and plots a view of the manifold into which individual points can be plotted later.`m.embed(v)`

- function handle, for given points`v`

on the manifold returns a points in the same (2D or 3D) space as`m.plot`

that are used to visualize the points`v`

.

The easiest way to find out the correct dimensions is to check out `manifold_moebius.m`

or `manifold_Rn.m`

. The solver will display an error if any of the dimensions are inconsistent with the `m.mdims`

struct.

If you have any questions or need help, feel free to send me an e-mail at Jan.Lellmann (at) web.de.

## Acknowledgments

- The computation of the area on the sphere uses parts of the SPHERE_DELAUNAY and SPHERE_CVT libraries written by John Burkardt.