The Riemannian Chambolle-Pock algorithm

The Riemannian Chambolle—Pock is a generalization of the Chambolle—Pock algorithm Chambolle and Pock [CP11] It is also known as primal-dual hybrid gradient (PDHG) or primal-dual proximal splitting (PDPS) algorithm.

In order to minimize over $p∈\mathcal M$ the cost function consisting of

\[F(p) + G(Λ(p)),\]

over $p∈\mathcal M$

where $F:\mathcal M → \overline{ℝ}$, $G:\mathcal N → \overline{ℝ}$, and $Λ:\mathcal M →\mathcal N$. If the manifolds $\mathcal M$ or $\mathcal N$ are not Hadamard, it has to be considered locally only, that is on geodesically convex sets $\mathcal C \subset \mathcal M$ and $\mathcal D \subset\mathcal N$ such that $Λ(\mathcal C) \subset \mathcal D$.

The algorithm is available in four variants: exact versus linearized (see variant) as well as with primal versus dual relaxation (see relax). For more details, see Bergmann, Herzog, Silva Louzeiro, Tenbrinck and Vidal-Núñez [BHS+21]. In the following, the case of the exact, primal relaxed Riemannian Chambolle—Pock algorithm is described.

Given base points $m∈\mathcal C$, $n=Λ(m)∈\mathcal D$, initial primal and dual values $p^{(0)} ∈\mathcal C$, $ξ_n^{(0)} ∈T_n^*\mathcal N$, and primal and dual step sizes $\sigma_0$, $\tau_0$, relaxation $\theta_0$, as well as acceleration $\gamma$.

As an initialization, perform $\bar p^{(0)} \gets p^{(0)}$.

The algorithm performs the steps $k=1,…,$ (until a StoppingCriterion is fulfilled)

  1. \[ξ^{(k+1)}_n = \operatorname{prox}_{\tau_k G_n^*}\Bigl(ξ_n^{(k)} + \tau_k \bigl(\log_n Λ (\bar p^{(k)})\bigr)^\flat\Bigr)\]

  2. \[p^{(k+1)} = \operatorname{prox}_{\sigma_k F}\biggl(\exp_{p^{(k)}}\Bigl( \operatorname{PT}_{p^{(k)}\gets m}\bigl(-\sigma_k DΛ(m)^*[ξ_n^{(k+1)}]\bigr)^\sharp\Bigr)\biggr)\]

  3. Update
    • $\theta_k = (1+2\gamma\sigma_k)^{-\frac{1}{2}}$
    • $\sigma_{k+1} = \sigma_k\theta_k$
    • $\tau_{k+1} = \frac{\tau_k}{\theta_k}$
  4. \[\bar p^{(k+1)} = \exp_{p^{(k+1)}}\bigl(-\theta_k \log_{p^{(k+1)}} p^{(k)}\bigr)\]

Furthermore you can exchange the exponential map, the logarithmic map, and the parallel transport by a retraction, an inverse retraction, and a vector transport.

Finally you can also update the base points $m$ and $n$ during the iterations. This introduces a few additional vector transports. The same holds for the case $Λ(m^{(k)})\neq n^{(k)}$ at some point. All these cases are covered in the algorithm.

Manopt.ChambollePockFunction
ChambollePock(M, N, f, p, X, m, n, prox_G, prox_G_dual, get_adjoint_linear_operator; kwargs...)
ChambollePock!(M, N, f, p, X, m, n, prox_G, prox_G_dual, get_adjoint_linear_operator; kwargs...)

Perform the Riemannian Chambolle-Pock algorithm.

Given a cost function $\mathcal{E}: \mathcal{M} → ℝ$ of the form

\[\mathcal{E}(p) = F(p) + G( Λ(p) ),\]

where $F:\mathcal{M} → ℝ$, $G:\mathcal{N} → ℝ$, and $Λ:\mathcal{M} → \mathcal{N}$.

This can be done inplace of $p$.

Input parameters

  • M::AbstractManifold: a Riemannian manifold $\mathcal{M}$
  • N::AbstractManifold: a Riemannian manifold $\mathcal{M}$
  • p::P: a point on the manifold $\mathcal{M}$
  • X::T: a tangent vector at the point $p$ on the manifold $\mathcal{M}$
  • m::P: a point on the manifold $\mathcal{M}$
  • n::P: a point on the manifold $\mathcal{N}$
  • adjoint_linearized_operator: the adjoint $DΛ^*$ of the linearized operator $DΛ: T_{m}\mathcal{M} → T_{Λ(m)}\mathcal{N}$
  • prox_F, prox_G_Dual: the proximal maps of $F$ and $G^\ast_n$

By default, this performs the exact Riemannian Chambolle Pock algorithm, see the optional parameter for their linearized variant.

For more details on the algorithm, see [BHS+21].

Keyword Arguments

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.

source
Manopt.ChambollePock!Function
ChambollePock(M, N, f, p, X, m, n, prox_G, prox_G_dual, get_adjoint_linear_operator; kwargs...)
ChambollePock!(M, N, f, p, X, m, n, prox_G, prox_G_dual, get_adjoint_linear_operator; kwargs...)

Perform the Riemannian Chambolle-Pock algorithm.

Given a cost function $\mathcal{E}: \mathcal{M} → ℝ$ of the form

\[\mathcal{E}(p) = F(p) + G( Λ(p) ),\]

where $F:\mathcal{M} → ℝ$, $G:\mathcal{N} → ℝ$, and $Λ:\mathcal{M} → \mathcal{N}$.

This can be done inplace of $p$.

Input parameters

  • M::AbstractManifold: a Riemannian manifold $\mathcal{M}$
  • N::AbstractManifold: a Riemannian manifold $\mathcal{M}$
  • p::P: a point on the manifold $\mathcal{M}$
  • X::T: a tangent vector at the point $p$ on the manifold $\mathcal{M}$
  • m::P: a point on the manifold $\mathcal{M}$
  • n::P: a point on the manifold $\mathcal{N}$
  • adjoint_linearized_operator: the adjoint $DΛ^*$ of the linearized operator $DΛ: T_{m}\mathcal{M} → T_{Λ(m)}\mathcal{N}$
  • prox_F, prox_G_Dual: the proximal maps of $F$ and $G^\ast_n$

By default, this performs the exact Riemannian Chambolle Pock algorithm, see the optional parameter for their linearized variant.

For more details on the algorithm, see [BHS+21].

Keyword Arguments

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.

source

Objective

Manopt.PrimalDualManifoldObjectiveType
PrimalDualManifoldObjective <: AbstractPrimalDualManifoldObjective

Describes an Objective linearized or exact Chambolle-Pock algorithm, cf. [BHS+21], [CP11]

Fields

  • cost: $F + G(Λ(⋅))$ to evaluate interim cost function values
  • linearized_forward_operator!: linearized operator for the forward operation in the algorithm $DΛ$
  • linearized_adjoint_operator!: the adjoint differential $(DΛ)^* : \mathcal{N} → T\mathcal{M}$
  • prox_f!: the proximal map belonging to $f$
  • prox_G_dual!: the proximal map belonging to $g_n^*$
  • Λ!: the forward operator (if given) $Λ: \mathcal{M} → \mathcal{N}$

Either the linearized operator $DΛ$ or $Λ$ are required usually.

Constructor

PrimalDualManifoldObjective(cost, prox_f, prox_G_dual, adjoint_linearized_operator;    linearized_forward_operator::Union{Function,Missing}=missing,    Λ::Union{Function,Missing}=missing,    evaluation::AbstractEvaluationType = AllocatingEvaluation(),)

Using the evaluation= keyword can be used to specify that all functions work in-place instead of the default allocating one.

source

State

Manopt.ChambollePockStateType
ChambollePockState <: AbstractPrimalDualSolverState

stores all options and variables within a linearized or exact Chambolle Pock.

Fields

  • callbacks::D: provided callback functions given as a dictionary with symbols as keys
  • acceleration::R: acceleration factor
  • dual_stepsize::R: proximal parameter of the dual prox
  • inverse_retraction_method::AbstractInverseRetractionMethod: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses
  • inverse_retraction_method_dual::AbstractInverseRetractionMethod: an inverse retraction $\operatorname{retr}^{-1}$ to use, see the section on retractions and their inverses
  • m::P: base point on $\mathcal{M}$
  • n::Q: base point on $\mathcal{N}$
  • p::P: an initial point on $p^{(0)} ∈ \mathcal{M}$
  • pbar::P: the relaxed iterate used in the next dual update step (when using :primal relaxation)
  • primal_stepsize::R: proximal parameter of the primal prox
  • X::T: an initial tangent vector $X^{(0)} ∈ T_{p^{(0)}}\mathcal{M}$
  • Xbar::T: the relaxed iterate used in the next primal update step (when using :dual relaxation)
  • relaxation::R: relaxation in the primal relaxation step (to compute pbar:
  • relax::Symbol: which variable to relax (:primalor:dual`:
  • retraction_method::AbstractRetractionMethod: a retraction $\operatorname{retr}$ to use, see the section on retractions
  • stop::StoppingCriterion: a functor indicating that the stopping criterion is fulfilled
  • variant: whether to perform an :exact or :linearized Chambolle-Pock
  • update_primal_base: function (pr, st, k) -> m to update the primal base
  • update_dual_base: function (pr, st, k) -> n to update the dual base
  • vector_transport_method::AbstractVectorTransportMethod: a vector transport $\mathcal T_{⋅←⋅}$ to use, see the section on vector transports
  • vector_transport_method_dual::AbstractVectorTransportMethod: a vector transport $\mathcal T_{⋅←⋅}$ to use, see the section on vector transports

Here, P is a point type on $\mathcal{M}$, T its tangent vector type, Q a point type on $\mathcal{N}$, and R<:Real is a real number type

where for the last two the functions a AbstractManoptProblemp, AbstractManoptSolverStateo and the current iterate i are the arguments. If you activate these to be different from the default identity, you have to provide p.Λ for the algorithm to work (which might be missing in the linearized case).

Constructor

ChambollePockState(M::AbstractManifold, N::AbstractManifold;    kwargs...) where {P, Q, T, R <: Real}

Keyword arguments

if Manifolds.jl is loaded, N is also a keyword argument and set to TangentBundle(M) by default.

source

Useful terms

Manopt.primal_residualFunction
primal_residual(tmp::TwoManifoldProblem, apds::AbstractPrimalDualSolverState, p_old, X_old, n_old)

Compute the primal residual at current iterate $k$ given the necessary values $p_{k-1}, X_{k-1}$, and $n_{k-1}$ from the previous iterate.

\[\lVert \frac{1}{σ}\operatorname{retr}^{-1}_{x_{k}}x_{k-1} - V_{x_k←m_k} \bigl( DΛ^*(m_k)\bigl[V_{n_k← n_{k-1}}X_{k-1} - X_k \bigr]\bigr) \rVert\]

where $V_{⋅←⋅}$ is the vector transport used in the ChambollePockState.

source
Manopt.dual_residualFunction
dual_residual(tmp::TwoManifoldProblem, apds::AbstractPrimalDualSolverState, p_old, X_old, n_old)

Compute the dual residual at current iterate $k$ given the necessary values $p_{k-1}, X_{k-1}$, and $n_{k-1}$ from the previous iterate. The formula is slightly different depending on the apds.variant used:

For the :linearized it reads

\[\lVert \frac{1}{τ}\bigl( V_{n_{k}← n_{k-1}}(X_{k-1}) - X_k \bigr ) - DΛ(m_k)\bigl[ V_{m_k← x_k}\operatorname{retr}^{-1}_{x_{k}}(x_{k-1})\bigr] \rVert\]

and for the :exact variant

\[\lVert \frac{1}{τ} V_{n_{k}← n_{k-1}}(X_{k-1}) - \operatorname{retr}^{-1}_{n_{k}}\bigl( Λ(\operatorname{retr}_{m_{k}}(V_{m_k← x_k}\operatorname{retr}^{-1}_{x_{k}}x_{k-1}))\bigr) \rVert\]

where in both cases $V_{⋅←⋅}$ is the vector transport used in the ChambollePockState.

source

Solver specific debug output

DebugDualBaseIterate, DebugDualBaseChange, DebugPrimalBaseIterate, DebugPrimalBaseChange, DebugDualChange, DebugDualIterate, DebugDualResidual, DebugPrimalChange, DebugPrimalIterate, DebugPrimalResidual DebugPrimalDualResidual

Solver specific recording actions

Internals

Manopt.update_prox_parameters!Function
update_prox_parameters!(o)

update the prox parameters as described in Algorithm 2 of [CP11],

  1. $θ_{n} = \frac{1}{\sqrt{1+2γτ_n}}$
  2. $τ_{n+1} = θ_nτ_n$
  3. $σ_{n+1} = \frac{σ_n}{θ_n}$
source

Technical details

The ChambollePock solver requires the following functions of a manifold to be available for both the manifold $\mathcal M$and $\mathcal N$

  • 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= or retraction_method_dual= (for $\mathcal N$) 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= or inverse_retraction_method_dual= (for $\mathcal N$) does not have to be specified.
  • A vector_transport_to!M, Y, p, X, q); it is recommended to set the default_vector_transport_method to a favourite retraction. If this default is set, a vector_transport_method= or vector_transport_method_dual= (for $\mathcal N$) does not have to be specified.
  • A copyto!(M, q, p) and copy(M,p) for points.

Literature

[BHS+21]
R. Bergmann, R. Herzog, M. Silva Louzeiro, D. Tenbrinck and J. Vidal-Núñez. Fenchel duality theory and a primal-dual algorithm on Riemannian manifolds. Foundations of Computational Mathematics 21, 1465–1504 (2021), arXiv:1908.02022.
[CP11]
A. Chambolle and T. Pock. A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision 40, 120–145 (2011).