The Riemannian Chambolle-Pock Algorithm

The Riemannian Chambolle–Pock is a generalization of the Chambolle–Pock algorithm[ChambollePock2011]. It is also known as primal dual hybrig gradient (PDHG) or primal dual proximal splitting (PDPS) algorithm.

In order to minimize over p\in\mathcal M§ the cost function consisting of

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

where $F:\mathcal M \to \overline{\mathbb R}$, $G:\mathcal N \to \overline{\mathbb R}$, and $\Lambda:\mathcal M \to\mathcal N$. If the manifolds $\mathcal M$ or $\mathcal N$ are not Hadamard, it has to be considered locally, i.e. on geodesically convex sets $\mathcal C \subset \mathcal M$ and $\mathcal D \subset\mathcal N$ such that $\Lambda(\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 [BergmannHerzogSilvaLouzeiroTenbrinckVidalNunez2020]. In the following we note the case of the exact, primal relaxed Riemannian Chambolle–Pock algorithm.

Given base points $m\in\mathcal C$, $n=\Lambda(m)\in\mathcal D$, initial primal and dual values $p^{(0)} \in \mathcal C$, $\xi_n^{(0)} \in 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 algorithms performs the steps $k=1,\ldots,$ (until a StoppingCriterion is fulfilled with)

  1. \[\xi^{(k+1)}_n = \operatorname{prox}_{\tau_k G_n^*}\Bigl(\xi_n^{(k)} + \tau_k \bigl(\log_n \Lambda (\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\Lambda(m)^*[\xi_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 in verse 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 that $\Lambda(m^{(k)})\neq n^{(k)}$ at some point. All these cases are covered in the algorithm.

Manopt.ChambollePockFunction
ChambollePock(M, N, cost, x0, ξ0, m, n, prox_F, prox_G_dual, forward_operator, adjoint_DΛ)

Perform the Riemannian Chambolle–Pock algorithm.

Given a cost function $\mathcal E\colon\mathcal M \to ℝ$ of the form

\[\mathcal E(x) = F(x) + G( Λ(x) ),\]

where $F\colon\mathcal M \to ℝ$, $G\colon\mathcal N \to ℝ$, and $\Lambda\colon\mathcal M \to \mathcal N$. The remaining input parameters are

  • x,ξ primal and dual start points $x\in\mathcal M$ and $\xi\in T_n\mathcal N$
  • m,n base points on $\mathcal M$ and $\mathcal N$, respectively.
  • forward_operator the operator $Λ(⋅)$ or its linearization $DΛ(⋅)[⋅]$, depending on whether :exact or :linearized is chosen.
  • adjoint_linearized_operator the adjoint $DΛ^*$ of the linearized operator $DΛ(m)\colon T_{m}\mathcal M \to 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 opional parameter for ther linearized variant.

For more details on the algorithm, see[BergmannHerzogSilvaLouzeiroTenbrinckVidalNunez2020].

Optional Parameters

  • acceleration – (0.05)
  • dual_stepsize – (1/sqrt(8)) proximnal parameter of the primal prox
  • Λ (missing) the exact operator, that is required if the forward operator is linearized; missing indicates, that the forward operator is exact.
  • primal_stepsize – (1/sqrt(8)) proximnal parameter of the dual prox
  • relaxation – (1.)
  • relax – (:primal) whether to relax the primal or dual
  • variant - (:exact if Λ is missing, otherwise :linearized) variant to use. Note that this changes the arguments the forward_operator will be called.
  • stopping_criterion – (stopAtIteration(100)) a StoppingCriterion
  • update_primal_base – (missing) function to update m (identity by default/missing)
  • update_dual_base – (missing) function to update n (identity by default/missing)
  • retraction_method – (ExponentialRetraction()) the rectraction to use
  • inverse_retraction_method - (LogarithmicInverseRetraction()) an inverse retraction to use.
  • vector_transport_method - (ParallelTransport()) a vector transport to use
source
Manopt.ChambollePock!Function
ChambollePock(M, N, cost, x, ξ, m, n, prox_F, prox_G_dual, forward_operator, adjoint_DΛ)

Perform the Riemannian Chambolle–Pock algorithm in place of x, ξ, and potenitally m, n if they are not fixed. See ChambollePock for details and optional parameters.

source

Problem & Options

Manopt.PrimalDualProblemType
PrimalDualProblem {mT <: Manifold, nT <: Manifold} <: PrimalDualProblem} <: Problem

Describes a Problem for the linearized Chambolle-Pock algorithm.

Fields

  • M, N – two manifolds $\mathcal M$, $\mathcal N$
  • cost $F + G(Λ(⋅))$ to evaluate interims cost function values
  • forward_oprator the operator for the forward operation in the algorthm, either $Λ$ (exact) or $DΛ$ (linearized).
  • linearized_adjoint_operator The adjoint differential $(DΛ)^* \colon \mathcal N \to T\mathcal M$
  • prox_F the proximal map belonging to $f$
  • prox_G_dual the proximal map belonging to $g_n^*$
  • Λ – (fordward_operator) for the linearized variant, this has to be set to the exact forward operator. This operator is required in several variants of the linearized algorithm. Since the exact variant is the default, Λ is by default set to forward_operator.

Constructor

LinearizedPrimalDualProblem(M, N, cost, prox_F, prox_G_dual, forward_operator, adjoint_linearized_operator,Λ=forward_operator)
source
Manopt.PrimalDualOptionsType
PrimalDualOptions

A general type for all primal dual based options to be used within primal dual based algorithms

source
Manopt.ChambollePockOptionsType
ChambollePockOptions <: PrimalDualOptions

stores all options and variables within a linearized or exact Chambolle Pock. The following list provides the order for the constructor, where the previous iterates are initialized automatically and values with a default may be left out.

  • m - base point on $ \mathcal M $
  • n - base point on $ \mathcal N $
  • x - an initial point on $x^{(0)} \in \mathcal M$ (and its previous iterate)
  • ξ - an initial tangent vector $\xi^{(0)}\in T^*\mathcal N$ (and its previous iterate)
  • xbar - the relaxed iterate used in the next dual update step (when using :primal relaxation)
  • ξbar - the relaxed iterate used in the next primal update step (when using :dual relaxation)
  • Θ – factor to damp the helping $\tilde x$
  • primal_stepsize – (1/sqrt(8)) proximal parameter of the primal prox
  • dual_stepsize – (1/sqrt(8)) proximnal parameter of the dual prox
  • acceleration – (0.) acceleration factor due to Chambolle & Pock
  • relaxation – (1.) relaxation in the primal relaxation step (to compute xbar)
  • relax – (_primal) which variable to relax (:primal or :dual)
  • stop - a StoppingCriterion
  • type – (exact) whether to perform an :exact or :linearized Chambolle-Pock
  • update_primal_base ((p,o,i) -> o.m) function to update the primal base
  • update_dual_base ((p,o,i) -> o.n) function to update the dual base
  • retraction_method – (ExponentialRetraction()) the rectraction to use
  • inverse_retraction_method - (LogarithmicInverseRetraction()) an inverse retraction to use.
  • vector_transport_method - (ParallelTransport()) a vector transport to use

where for the last two the functions a Problemp, Optionso 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

ChambollePockOptions(m::P, n::Q, x::P, ξ::T, primal_stepsize::Float64, dual_stepsize::Float64;
    acceleration::Float64 = 0.0,
    relaxation::Float64 = 1.0,
    relax::Symbol = :primal,
    stopping_criterion::StoppingCriterion = StopAfterIteration(300),
    variant::Symbol = :exact,
    update_primal_base::Union{Function,Missing} = missing,
    update_dual_base::Union{Function,Missing} = missing,
    retraction_method = ExponentialRetraction(),
    inverse_retraction_method = LogarithmicInverseRetraction(),
    vector_transport_method = ParallelTransport(),
)
source

Useful Terms

Manopt.primal_residualFunction
primal_residual(p, o, x_old, ξ_old, n_old)

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

\[\Bigl\lVert \frac{1}{σ}\operatorname{retr}^{-1}_{x_{k}}x_{k-1} - V_{x_k\gets m_k}\bigl(DΛ^*(m_k)\bigl[V_{n_k\gets n_{k-1}}ξ_{k-1} - ξ_k \bigr] \Bigr\rVert\]

where $V_{\cdot\gets\cdot}$ is the vector transport used in the ChambollePockOptions

source
Manopt.dual_residualFunction
dual_residual(p, o, x_old, ξ_old, n_old)

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

For the :lineaized it reads

\[\Bigl\lVert \frac{1}{τ}\bigl( V_{n_{k}\gets n_{k-1}}(ξ_{k-1}) - ξ_k \bigr) - DΛ(m_k)\bigl[ V_{m_k\gets x_k}\operatorname{retr}^{-1}_{x_{k}}x_{k-1} \bigr] \Bigr\rVert\]

and for the :exact variant

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

where in both cases $V_{\cdot\gets\cdot}$ is the vector transport used in the ChambollePockOptions.

source

Debug

Manopt.DebugDualChangeType
DebugDualChange(opts...)

Print the change of the dual variable, similar to DebugChange, see their constructors for detail, but with a different calculation of the change, since the dual variable lives in (possibly different) tangent spaces.

source
Manopt.DebugDualResidualType
DebugDualResidual <: DebugAction

A Debug action to print the dual residual. The constructor accepts a printing function and some (shared) storage, which should at least record :x, and :n.

source
Manopt.DebugPrimalResidualType
DebugPrimalResidual <: DebugAction

A Debug action to print the primal residual. The constructor accepts a printing function and some (shared) storage, which should at least record :x, and :n.

source
Manopt.DebugPrimalDualResidualType
DebugPrimalDualResidual <: DebugAction

A Debug action to print the primaldual residual. The constructor accepts a printing function and some (shared) storage, which should at least record :x, and :n.

source

Record

Manopt.RecordDualChangeFunction
RecordDualChange()

Create the action either with a given (shared) Storage, which can be set to the values Tuple, if that is provided).

source

Internals

Manopt.update_prox_parameters!Function
update_prox_parameters!(o)

update the prox parameters as described in Algorithm 2 of Chambolle, Pock, 2010, i.e.

  1. $\theta_{n} = \frac{1}{\sqrt{1+2\gamma\tau_n}}$
  2. $\tau_{n+1} = \theta_n\tau_n$
  3. $\sigma_{n+1} = \frac{\sigma_n}{\theta_n}$
source
  • BergmannHerzogSilvaLouzeiroTenbrinckVidalNunez2020

    R. Bergmann, R. Herzog, M. Silva Louzeiro, D. Tenbrinck, J. Vidal-Núñez: Fenchel Duality Theory and a Primal-Dual Algorithm on Riemannian Manifolds, arXiv: 1908.02022 accepted for publication in Foundations of Computational Mathematics

  • ChambollePock2011

    A. Chambolle, T. Pock: A first-order primal-dual algorithm for convex problems with applications to imaging, Journal of Mathematical Imaging and Vision 40(1), 120–145, 2011. doi: 10.1007/s10851-010-0251-1