# Douglas—Rachford algorithm

The (Parallel) Douglas—Rachford ((P)DR) algorithm was generalized to Hadamard manifolds in [BPS16].

The aim is to minimize the sum

\[f(p) = g(p) + h(p)\]

on a manifold, where the two summands have proximal maps $\operatorname{prox}_{λ g}, \operatorname{prox}_{λ h}$ that are easy to evaluate (maybe in closed form, or not too costly to approximate). Further, define the reflection operator at the proximal map as

\[\operatorname{refl}_{λ g}(p) = \operatorname{retr}_{\operatorname{prox}_{λ g}(p)} \bigl( -\operatorname{retr}^{-1}_{\operatorname{prox}_{λ g}(p)} p \bigr).\]

Let $\alpha_k ∈ [0,1]$ with $\sum_{k ∈ ℕ} \alpha_k(1-\alpha_k) = \infty$ and $λ > 0$ (which might depend on iteration $k$ as well) be given.

Then the (P)DRA algorithm for initial data $p^{(0)} ∈ \mathcal M$ as

## Initialization

Initialize $q^{(0)} = p^{(0)}$ and $k=0$

## Iteration

Repeat until a convergence criterion is reached

- Compute $r^{(k)} = \operatorname{refl}_{λ g}\operatorname{refl}_{λ h}(q^{(k)})$
- Within that operation, store $p^{(k+1)} = \operatorname{prox}_{λ h}(q^{(k)})$ which is the prox the inner reflection reflects at.
- Compute $q^{(k+1)} = g(\alpha_k; q^{(k)}, r^{(k)})$, where $g$ is a curve approximating the shortest geodesic, provided by a retraction and its inverse
- Set $k = k+1$

## Result

The result is given by the last computed $p^{(K)}$ at the last iterate $K$.

For the parallel version, the first proximal map is a vectorial version where in each component one prox is applied to the corresponding copy of $t_k$ and the second proximal map corresponds to the indicator function of the set, where all copies are equal (in $\mathcal M^n$, where $n$ is the number of copies), leading to the second prox being the Riemannian mean.

## Interface

`Manopt.DouglasRachford`

— Function```
DouglasRachford(M, f, proxes_f, p)
DouglasRachford(M, mpo, p)
DouglasRachford!(M, f, proxes_f, p)
DouglasRachford!(M, mpo, p)
```

Compute the Douglas-Rachford algorithm on the manifold $\mathcal M$, starting from `p`

`given the (two) proximal maps`

proxes_f`, see [BPS16].

For $k>2$ proximal maps, the problem is reformulated using the parallel Douglas Rachford: a vectorial proximal map on the power manifold $\mathcal M^k$ is introduced as the first proximal map and the second proximal map of the is set to the `mean`

(Riemannian center of mass). This hence also boils down to two proximal maps, though each evaluates proximal maps in parallel, that is, component wise in a vector.

The parallel Douglas Rachford does not work in-place for now, since while creating the new staring point `p'`

on the power manifold, a copy of `p`

Is created

If you provide a `ManifoldProximalMapObjective`

`mpo`

instead, the proximal maps are kept unchanged.

**Input**

`M::`

`AbstractManifold`

`: a Riemannian manifold $\mathcal M$`

`f`

: a cost function $f: \mathcal M→ ℝ$ implemented as`(M, p) -> v`

`proxes_f`

: functions of the form`(M, λ, p)-> q`

performing a proximal maps, where`λ`

denotes the proximal parameter, for each of the summands of`F`

. These can also be given in the`InplaceEvaluation`

variants`(M, q, λ p) -> q`

computing in place of`q`

.`p`

: a point on the manifold $\mathcal M$

**Keyword arguments**

`α= k -> 0.9`

: relaxation of the step from old to new iterate, to be precise $p^{(k+1)} = g(α_k; p^{(k)}, q^{(k)})$, where $q^{(k)}$ is the result of the double reflection involved in the DR algorithm and $g$ is a curve induced by the retraction and its inverse.`evaluation=`

`AllocatingEvaluation`

`()`

: specify whether the functions that return an array, for example a point or a tangent vector, work by allocating its result (`AllocatingEvaluation`

) or whether they modify their input argument to return the result therein (`InplaceEvaluation`

). Since usually the first argument is the manifold, the modified argument is the second.`inverse_retraction_method=`

`default_inverse_retraction_method`

`(M, typeof(p))`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses This is used both in the relaxation step as well as in the reflection, unless you set`R`

yourself.`λ= k -> 1.0`

: function to provide the value for the proximal parameter $λ_k$`R=reflect(!)`

: method employed in the iteration to perform the reflection of`p`

at the prox of`p`

. This uses by default`reflect`

or`reflect!`

depending on`reflection_evaluation`

and the retraction and inverse retraction specified by`retraction_method`

and`inverse_retraction_method`

, respectively.`reflection_evaluation`

: (`AllocatingEvaluation`

whether`R`

works in-place or allocating`retraction_method=`

`default_retraction_method`

`(M, typeof(p))`

: a retraction $\operatorname{retr}$ to use, see the section on retractions This is used both in the relaxation step as well as in the reflection, unless you set`R`

yourself.`stopping_criterion=`

`StopAfterIteration`

`(200)`

`|`

`StopWhenChangeLess`

`(1e-5)`

: a functor indicating that the stopping criterion is fulfilled`parallel=false`

: indicate whether to use a parallel Douglas-Rachford or not.

All other keyword arguments are passed to `decorate_state!`

for state decorators or `decorate_objective!`

for objective decorators, respectively.

**Output**

The obtained approximate minimizer $p^*$. To obtain the whole final state of the solver, see `get_solver_return`

for details, especially the `return_state=`

keyword.

`DouglasRachford(M, f, proxes_f, p; kwargs...)`

a doc string with some math $t_{k+1} = g(α_k; t_k, s_k)$

`Manopt.DouglasRachford!`

— Function```
DouglasRachford(M, f, proxes_f, p)
DouglasRachford(M, mpo, p)
DouglasRachford!(M, f, proxes_f, p)
DouglasRachford!(M, mpo, p)
```

Compute the Douglas-Rachford algorithm on the manifold $\mathcal M$, starting from `p`

`given the (two) proximal maps`

proxes_f`, see [BPS16].

For $k>2$ proximal maps, the problem is reformulated using the parallel Douglas Rachford: a vectorial proximal map on the power manifold $\mathcal M^k$ is introduced as the first proximal map and the second proximal map of the is set to the `mean`

(Riemannian center of mass). This hence also boils down to two proximal maps, though each evaluates proximal maps in parallel, that is, component wise in a vector.

The parallel Douglas Rachford does not work in-place for now, since while creating the new staring point `p'`

on the power manifold, a copy of `p`

Is created

If you provide a `ManifoldProximalMapObjective`

`mpo`

instead, the proximal maps are kept unchanged.

**Input**

`M::`

`AbstractManifold`

`: a Riemannian manifold $\mathcal M$`

`f`

: a cost function $f: \mathcal M→ ℝ$ implemented as`(M, p) -> v`

`proxes_f`

: functions of the form`(M, λ, p)-> q`

performing a proximal maps, where`λ`

denotes the proximal parameter, for each of the summands of`F`

. These can also be given in the`InplaceEvaluation`

variants`(M, q, λ p) -> q`

computing in place of`q`

.`p`

: a point on the manifold $\mathcal M$

**Keyword arguments**

`α= k -> 0.9`

: relaxation of the step from old to new iterate, to be precise $p^{(k+1)} = g(α_k; p^{(k)}, q^{(k)})$, where $q^{(k)}$ is the result of the double reflection involved in the DR algorithm and $g$ is a curve induced by the retraction and its inverse.`evaluation=`

`AllocatingEvaluation`

`()`

: specify whether the functions that return an array, for example a point or a tangent vector, work by allocating its result (`AllocatingEvaluation`

) or whether they modify their input argument to return the result therein (`InplaceEvaluation`

). Since usually the first argument is the manifold, the modified argument is the second.`inverse_retraction_method=`

`default_inverse_retraction_method`

`(M, typeof(p))`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses This is used both in the relaxation step as well as in the reflection, unless you set`R`

yourself.`λ= k -> 1.0`

: function to provide the value for the proximal parameter $λ_k$`R=reflect(!)`

: method employed in the iteration to perform the reflection of`p`

at the prox of`p`

. This uses by default`reflect`

or`reflect!`

depending on`reflection_evaluation`

and the retraction and inverse retraction specified by`retraction_method`

and`inverse_retraction_method`

, respectively.`reflection_evaluation`

: (`AllocatingEvaluation`

whether`R`

works in-place or allocating`retraction_method=`

`default_retraction_method`

`(M, typeof(p))`

: a retraction $\operatorname{retr}$ to use, see the section on retractions This is used both in the relaxation step as well as in the reflection, unless you set`R`

yourself.`stopping_criterion=`

`StopAfterIteration`

`(200)`

`|`

`StopWhenChangeLess`

`(1e-5)`

: a functor indicating that the stopping criterion is fulfilled`parallel=false`

: indicate whether to use a parallel Douglas-Rachford or not.

All other keyword arguments are passed to `decorate_state!`

for state decorators or `decorate_objective!`

for objective decorators, respectively.

**Output**

The obtained approximate minimizer $p^*$. To obtain the whole final state of the solver, see `get_solver_return`

for details, especially the `return_state=`

keyword.

## State

`Manopt.DouglasRachfordState`

— Type`DouglasRachfordState <: AbstractManoptSolverState`

Store all options required for the DouglasRachford algorithm,

**Fields**

`α`

: relaxation of the step from old to new iterate, to be precise $x^{(k+1)} = g(α(k); x^{(k)}, t^{(k)})$, where $t^{(k)}$ is the result of the double reflection involved in the DR algorithm`inverse_retraction_method::AbstractInverseRetractionMethod`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses`λ`

: function to provide the value for the proximal parameter during the calls`parallel`

: indicate whether to use a parallel Douglas-Rachford or not.`R`

: method employed in the iteration to perform the reflection of`x`

at the prox`p`

.`p::P`

: a point on the manifold $\mathcal M$storing the current iterate For the parallel Douglas-Rachford, this is not a value from the`PowerManifold`

manifold but the mean.`reflection_evaluation`

: whether`R`

works in-place or allocating`retraction_method::AbstractRetractionMethod`

: a retraction $\operatorname{retr}$ to use, see the section on retractions`s`

: the last result of the double reflection at the proximal maps relaxed by`α`

.`stop::StoppingCriterion`

: a functor indicating that the stopping criterion is fulfilled

**Constructor**

`DouglasRachfordState(M::AbstractManifold; kwargs...)`

**Input**

`M::`

`AbstractManifold`

`: a Riemannian manifold $\mathcal M$`

**Keyword arguments**

`α= k -> 0.9`

: relaxation of the step from old to new iterate, to be precise $x^{(k+1)} = g(α(k); x^{(k)}, t^{(k)})$, where $t^{(k)}$ is the result of the double reflection involved in the DR algorithm`inverse_retraction_method=`

`default_inverse_retraction_method`

`(M, typeof(p))`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses`λ= k -> 1.0`

: function to provide the value for the proximal parameter during the calls`p=`

`rand`

`(M)`

: a point on the manifold $\mathcal M$to specify the initial value`R=`

`reflect`

`(!)`

: method employed in the iteration to perform the reflection of`p`

at the prox of`p`

, which function is used depends on`reflection_evaluation`

.`reflection_evaluation=`

`AllocatingEvaluation`

`()`

) specify whether the reflection works in-place or allocating (default)`retraction_method=`

`default_retraction_method`

`(M, typeof(p))`

: a retraction $\operatorname{retr}$ to use, see the section on retractions`stopping_criterion=`

`StopAfterIteration`

`(300)`

: a functor indicating that the stopping criterion is fulfilled`parallel=false`

: indicate whether to use a parallel Douglas-Rachford or not.

For specific `DebugAction`

s and `RecordAction`

s see also Cyclic Proximal Point.

Furthermore, this solver has a short hand notation for the involved `reflect`

ion.

`Manopt.reflect`

— Function```
reflect(M, f, x; kwargs...)
reflect!(M, q, f, x; kwargs...)
```

reflect the point `x`

from the manifold `M`

at the point `f(x)`

of the function $f: \mathcal M → \mathcal M$, given by

\[ \operatorname{refl}_f(x) = \operatorname{refl}_{f(x)}(x),\]

Compute the result in `q`

.

see also `reflect`

`(M,p,x)`

, to which the keywords are also passed to.

```
reflect(M, p, x, kwargs...)
reflect!(M, q, p, x, kwargs...)
```

Reflect the point `x`

from the manifold `M`

at point `p`

, given by

\[\operatorname{refl}\]

where $\operatorname{retr}$ and $\operatorname{retr}^{-1}$ denote a retraction and an inverse retraction, respectively. This can also be done in place of `q`

.

**Keyword arguments**

`retraction_method=`

`default_retraction_method`

`(M, typeof(p))`

: a retraction $\operatorname{retr}$ to use, see the section on retractions`inverse_retraction_method=`

`default_inverse_retraction_method`

`(M, typeof(p))`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses

and for the `reflect!`

additionally

`X=`

`zero_vector`

`(M, p)`

: a tangent vector at the point $p$ on the manifold $\mathcal M$ as temporary memory to compute the inverse retraction in place. otherwise this is the memory that would be allocated anyways.

```
reflect(M, f, x; kwargs...)
reflect!(M, q, f, x; kwargs...)
```

reflect the point `x`

from the manifold `M`

at the point `f(x)`

of the function $f: \mathcal M → \mathcal M$, given by

\[ \operatorname{refl}_f(x) = \operatorname{refl}_{f(x)}(x),\]

Compute the result in `q`

.

see also `reflect`

`(M,p,x)`

, to which the keywords are also passed to.

```
reflect(M, p, x, kwargs...)
reflect!(M, q, p, x, kwargs...)
```

Reflect the point `x`

from the manifold `M`

at point `p`

, given by

\[\operatorname{refl}_p(q) = \operatorname{retr}_p(-\operatorname{retr}^{-1}_p q),\]

where $\operatorname{retr}$ and $\operatorname{retr}^{-1}$ denote a retraction and an inverse retraction, respectively.

This can also be done in place of `q`

.

**Keyword arguments**

`retraction_method=`

`default_retraction_method`

`(M, typeof(p))`

: a retraction $\operatorname{retr}$ to use, see the section on retractions`inverse_retraction_method=`

`default_inverse_retraction_method`

`(M, typeof(p))`

: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses

and for the `reflect!`

additionally

`X=zero_vector(M,p)`

: a temporary memory to compute the inverse retraction in place. otherwise this is the memory that would be allocated anyways.

## Technical details

The `DouglasRachford`

solver requires the following functions of a manifold to be available

- A
`retract!`

`(M, q, p, X)`

; it is recommended to set the`default_retraction_method`

to a favourite retraction. If this default is set, a`retraction_method=`

does not have to be specified. - An
`inverse_retract!`

`(M, X, p, q)`

; it is recommended to set the`default_inverse_retraction_method`

to a favourite retraction. If this default is set, a`inverse_retraction_method=`

does not have to be specified. - A `copyto!
`(M, q, p)`

and`copy`

`(M,p)`

for points.

By default, one of the stopping criteria is `StopWhenChangeLess`

, which requires

- An
`inverse_retract!`

`(M, X, p, q)`

; it is recommended to set the`default_inverse_retraction_method`

to a favourite retraction. If this default is set, a`inverse_retraction_method=`

or`inverse_retraction_method_dual=`

(for $\mathcal N$) does not have to be specified or the`distance`

`(M, p, q)`

for said default inverse retraction.