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) -> vproxes_f: functions of the form(M, λ, p)-> qperforming a proximal maps, whereλdenotes the proximal parameter, for each of the summands ofF. These can also be given in theInplaceEvaluationvariants(M, q, λ p) -> qcomputing in place ofq.p::P: a point on the manifold $\mathcal{M}$
Keyword arguments
callbacks::D = Dict{Symbol,Function}(): provided callback functions . A single function(symbol, problem, state, k)called in every hook, a (vector of) pairs:hook => functionor will be processed byprocess_callbacks_arg. As key you can either pass single symbol or an array of symbols to indicate a callback should be added in multiple placesα= 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::AbstractEvaluationType=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::AbstractInverseRetractionMethod=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 setRyourself.λ= k -> 1.0: function to provide the value for the proximal parameter $λ_k$R=reflect!: method employed in the iteration to perform the reflection ofpat the prox ofp. This uses by defaultreflectorreflect!depending onreflection_evaluationand the retraction and inverse retraction specified byretraction_methodandinverse_retraction_method, respectively.reflection_evaluation: (AllocatingEvaluationwhetherRworks in-place or allocatingretraction_method::AbstractRetractionMethod=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 setRyourself.stopping_criterion::StoppingCriterion=StopAfterIteration(200)|StopWhenChangeLess(1e-5): a functor indicating that the stopping criterion is fulfilledparallel=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.
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) -> vproxes_f: functions of the form(M, λ, p)-> qperforming a proximal maps, whereλdenotes the proximal parameter, for each of the summands ofF. These can also be given in theInplaceEvaluationvariants(M, q, λ p) -> qcomputing in place ofq.p::P: a point on the manifold $\mathcal{M}$
Keyword arguments
callbacks::D = Dict{Symbol,Function}(): provided callback functions . A single function(symbol, problem, state, k)called in every hook, a (vector of) pairs:hook => functionor will be processed byprocess_callbacks_arg. As key you can either pass single symbol or an array of symbols to indicate a callback should be added in multiple placesα= 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::AbstractEvaluationType=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::AbstractInverseRetractionMethod=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 setRyourself.λ= k -> 1.0: function to provide the value for the proximal parameter $λ_k$R=reflect!: method employed in the iteration to perform the reflection ofpat the prox ofp. This uses by defaultreflectorreflect!depending onreflection_evaluationand the retraction and inverse retraction specified byretraction_methodandinverse_retraction_method, respectively.reflection_evaluation: (AllocatingEvaluationwhetherRworks in-place or allocatingretraction_method::AbstractRetractionMethod=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 setRyourself.stopping_criterion::StoppingCriterion=StopAfterIteration(200)|StopWhenChangeLess(1e-5): a functor indicating that the stopping criterion is fulfilledparallel=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 <: AbstractManoptSolverStateStore 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 algorithmcallbacks::D: provided callback functions given as a dictionary with symbols as keysinverse_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 callsparallel: indicate whether to use a parallel Douglas-Rachford or not.R!: method employed in the iteration to perform the reflection ofxat the proxp.p::P: a point on the manifold $\mathcal{M}$ storing the current iterate For the parallel Douglas-Rachford, this is not a value from thePowerManifoldmanifold but the mean.reflection!: whetherRworks in-place or allocatingretraction_method::AbstractRetractionMethod: a retraction $\operatorname{retr}$ to use, see the section on retractionss: 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 algorithmcallbacks::D = Dict{Symbol,Function}(): provided callback functions given as a dictionary with symbols as keysinverse_retraction_method::AbstractInverseRetractionMethod=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 callsp::P =rand(M): a point on the manifold $\mathcal{M}$ to specify the initial valueR!=reflect: method employed in the iteration to perform the reflection ofpat the prox ofp, which always works in-place.reflection_evaluation=InplaceEvaluation()) specify whether the reflection works in-place (default) or allocatingretraction_method::AbstractRetractionMethod=default_retraction_method(M, typeof(p)): a retraction $\operatorname{retr}$ to use, see the section on retractionsstopping_criterion::StoppingCriterion=StopAfterIteration(300): a functor indicating that the stopping criterion is fulfilledparallel=false: indicate whether to use a parallel Douglas-Rachford or not.
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 thedefault_retraction_methodto a favourite retraction. If this default is set, aretraction_method=does not have to be specified. - An
inverse_retract!(M, X, p, q); it is recommended to set thedefault_inverse_retraction_methodto a favourite retraction. If this default is set, ainverse_retraction_method=does not have to be specified. - A
copyto!(M, q, p)andcopy(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 thedefault_inverse_retraction_methodto a favourite retraction. If this default is set, ainverse_retraction_method=orinverse_retraction_method_dual=(for $\mathcal N$) does not have to be specified or thedistance(M, p, q)for said default inverse retraction.
Literature
- [BPS16]
- R. Bergmann, J. Persch and G. Steidl. A parallel Douglas Rachford algorithm for minimizing ROF-like functionals on images with values in symmetric Hadamard manifolds. SIAM Journal on Imaging Sciences 9, 901–937 (2016), arXiv:1512.02814.