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
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)
- \[\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)\]
- \[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)\]
- 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}$
- \[\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.ChambollePock — FunctionChambollePock(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
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,nbase points on $\mathcal M$ and $\mathcal N$, respectively.forward_operatorthe operator $Λ(⋅)$ or its linearization $DΛ(⋅)[⋅]$, depending on whether:exactor:linearizedis chosen.adjoint_linearized_operatorthe adjoint $DΛ^*$ of the linearized operator $DΛ(m)\colon T_{m}\mathcal M \to T_{Λ(m)}\mathcal N$prox_F, prox_G_Dualthe proximal maps of $F$ and $G^\ast_n$
By default, this performs the exact Riemannian Chambolle Pock algorithm, see the opional parameter DΛ 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;missingindicates, that the forward operator is exact.primal_stepsize– (1/sqrt(8)) proximnal parameter of the dual proxrelaxation– (1.)relax– (:primal) whether to relax the primal or dualvariant- (:exactifΛis missing, otherwise:linearized) variant to use. Note that this changes the arguments theforward_operatorwill be called.stopping_criterion– (stopAtIteration(100)) aStoppingCriterionupdate_primal_base– (missing) function to updatem(identity by default/missing)update_dual_base– (missing) function to updaten(identity by default/missing)retraction_method– (ExponentialRetraction()) the rectraction to useinverse_retraction_method- (LogarithmicInverseRetraction()) an inverse retraction to use.vector_transport_method- (ParallelTransport()) a vector transport to use
Manopt.ChambollePock! — FunctionChambollePock(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.
Problem & Options
Manopt.PrimalDualProblem — TypePrimalDualProblem {mT <: Manifold, nT <: Manifold} <: PrimalDualProblem} <: ProblemDescribes 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 valuesforward_opratorthe operator for the forward operation in the algorthm, either $Λ$ (exact) or $DΛ$ (linearized).linearized_adjoint_operatorThe adjoint differential $(DΛ)^* \colon \mathcal N \to T\mathcal M$prox_Fthe proximal map belonging to $f$prox_G_dualthe 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 toforward_operator.
Constructor
LinearizedPrimalDualProblem(M, N, cost, prox_F, prox_G_dual, forward_operator, adjoint_linearized_operator,Λ=forward_operator)Manopt.PrimalDualOptions — TypePrimalDualOptionsA general type for all primal dual based options to be used within primal dual based algorithms
Manopt.ChambollePockOptions — TypeChambollePockOptions <: PrimalDualOptionsstores 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:primalrelaxation)ξbar- the relaxed iterate used in the next primal update step (when using:dualrelaxation)Θ– factor to damp the helping $\tilde x$primal_stepsize– (1/sqrt(8)) proximal parameter of the primal proxdual_stepsize– (1/sqrt(8)) proximnal parameter of the dual proxacceleration– (0.) acceleration factor due to Chambolle & Pockrelaxation– (1.) relaxation in the primal relaxation step (to computexbar)relax– (_primal) which variable to relax (:primalor:dual)stop- aStoppingCriteriontype– (exact) whether to perform an:exactor:linearizedChambolle-Pockupdate_primal_base((p,o,i) -> o.m) function to update the primal baseupdate_dual_base((p,o,i) -> o.n) function to update the dual baseretraction_method– (ExponentialRetraction()) the rectraction to useinverse_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(),
)Useful Terms
Manopt.primal_residual — Functionprimal_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.
where $V_{\cdot\gets\cdot}$ is the vector transport used in the ChambollePockOptions
Manopt.dual_residual — Functiondual_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
and for the :exact variant
where in both cases $V_{\cdot\gets\cdot}$ is the vector transport used in the ChambollePockOptions.
Debug
Manopt.DebugDualBaseIterate — FunctionDebugDualBaseIterate(io::IO=stdout)Print the dual base variable by using DebugEntry, see their constructors for detail. This method is further set display o.n.
Manopt.DebugDualBaseChange — FunctionDebugDualChange(a=StoreOptionsAction((:ξ)),io::IO=stdout)Print the change of the dual base variable by using DebugEntryChange, see their constructors for detail, on o.n.
Manopt.DebugPrimalBaseIterate — FunctionDebugPrimalBaseIterate(io::IO=stdout)Print the primal base variable by using DebugEntry, see their constructors for detail. This method is further set display o.m.
Manopt.DebugPrimalBaseChange — FunctionDebugPrimalBaseChange(a::StoreOptionsAction=StoreOptionsAction((:m)),io::IO=stdout)Print the change of the primal base variable by using DebugEntryChange, see their constructors for detail, on o.n.
Manopt.DebugDualChange — TypeDebugDualChange(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.
Manopt.DebugDualIterate — FunctionDebugDualIterate(e)Print the dual variable by using DebugEntry, see their constructors for detail. This method is further set display o.ξ.
Manopt.DebugDualResidual — TypeDebugDualResidual <: DebugActionA 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.
Manopt.DebugPrimalChange — FunctionDebugPrimalChange(opts...)Print the change of the primal variable by using DebugChange, see their constructors for detail.
Manopt.DebugPrimalIterate — FunctionDebugPrimalIterate(opts...)Print the change of the primal variable by using DebugIterate, see their constructors for detail.
Manopt.DebugPrimalResidual — TypeDebugPrimalResidual <: DebugActionA 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.
Manopt.DebugPrimalDualResidual — TypeDebugPrimalDualResidual <: DebugActionA 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.
Record
Manopt.RecordDualBaseIterate — FunctionRecordDualBaseIterate(n)Create an RecordAction that records the dual base point, i.e. RecordEntry of o.n.
Manopt.RecordDualBaseChange — FunctionRecordDualBaseChange(e)Create an RecordAction that records the dual base point change, i.e. RecordEntryChange of o.n with distance to the last value to store a value.
Manopt.RecordDualChange — FunctionRecordDualChange()Create the action either with a given (shared) Storage, which can be set to the values Tuple, if that is provided).
Manopt.RecordDualIterate — FunctionRecordDualIterate(ξ)Create an RecordAction that records the dual base point, i.e. RecordEntry of o.ξ, so .
Manopt.RecordPrimalBaseIterate — FunctionRecordPrimalBaseIterate(x)Create an RecordAction that records the primal base point, i.e. RecordEntry of o.m.
Manopt.RecordPrimalBaseChange — FunctionRecordPrimalBaseChange()Create an RecordAction that records the primal base point change, i.e. RecordEntryChange of o.m with distance to the last value to store a value.
Manopt.RecordPrimalChange — FunctionRecordPrimalChange(a)Create an RecordAction that records the primal value change, i.e. RecordChange, since we just redord the change of o.x.
Manopt.RecordPrimalIterate — FunctionRecordDualBaseIterate(x)Create an RecordAction that records the dual base point, i.e. RecordIterate, i.e. o.x.
Internals
Manopt.update_prox_parameters! — Functionupdate_prox_parameters!(o)update the prox parameters as described in Algorithm 2 of Chambolle, Pock, 2010, i.e.
- $\theta_{n} = \frac{1}{\sqrt{1+2\gamma\tau_n}}$
- $\tau_{n+1} = \theta_n\tau_n$
- $\sigma_{n+1} = \frac{\sigma_n}{\theta_n}$
- 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