---
title: Riemannian Block Coordinate Descent
url: https://www.emergentmind.com/topics/riemannian-block-coordinate-descent-rbcd
type: topic
---

# Riemannian Block Coordinate Descent

Riemannian Block Coordinate Descent (RBCD) denotes a family of block-wise optimization methods on Riemannian manifolds in which one updates one block variable, one tangent-space coordinate or block, or one tangent subspace at a time while maintaining feasibility through the exponential map or a retraction. In the recent literature, the term covers at least three closely related formulations: block coordinate descent on product manifolds \(M=M_1\times\cdots\times M_b\), tangent-subspace methods that replace Euclidean coordinate blocks by subspaces \(S_k\subseteq T_xM\), and coordinate descent algorithms on matrix manifolds built from manifold-specific tangent bases. The same acronym also names a specific fast, scalable algorithm for entropy-regularized Projection Robust Wasserstein distance on the Stiefel manifold [2305.14744] [1912.10627] [2406.02225] [2012.05199].

## 1. Problem formulations and basic definitions

A standard product-manifold formulation considers
\[
\min_{x\in M} F(x), \qquad M:=M_1\times\cdots\times M_b,
\]
where each block \(x_i\in M_i\subset \mathbb{R}^{n_i}\) and each \(M_i\) is a smooth embedded submanifold in the convergence-rate setting. The Euclidean partial gradient \(\nabla_iF(x)\) is projected onto \(T_{x_i}M_i\) to obtain the Riemannian partial gradient \(\widetilde{\nabla}_iF(x)\), and the full Riemannian gradient is \(\widetilde{\nabla}F(x)=[\widetilde{\nabla}_1F(x);\dots;\widetilde{\nabla}_bF(x)]\). For general closed sets, stationarity is expressed by the tangent cone condition
\[
\langle \nabla F(x),v\rangle \ge 0 \quad \forall v\in T_M(x),
\]
which reduces to \(\widetilde{\nabla}F(x)=0\) on smooth manifolds [2305.14744].

A second formulation, introduced as Tangent Subspace Descent (TSD), treats Euclidean coordinate blocks as tangent subspaces of a manifold. At iterate \(x_k\), one selects a tangent subspace \(S_k\subseteq T_{x_k}M\) with orthogonal projection \(P_{S_k}\), and updates by
\[
x_{k+1}=\Exp_{x_k}\big(-\alpha_k P_{S_k}\nabla f(x_k)\big),
\]
or, more generally, with a retraction \(R\),
\[
x_{k+1}=R_{x_k}\big(-\alpha_k P_{S_k}\nabla f(x_k)\big).
\]
In this formulation, RBCD is not tied to a product decomposition of variables; the “blocks” are tangent subspaces selected at the current iterate [1912.10627].

A third formulation specializes to matrix manifolds. There, coordinates are chosen with respect to a basis \(\{B_\ell\}\) spanning \(T_XM\), not necessarily orthonormal, and a block is a structured subset of entries or a tangent subspace. The generic update is
\[
x_{k+1}=R_{x_k}(\alpha_k\xi_k), \qquad \xi_k=P_{B_k}\,\mathrm{grad}\,f(x_k),
\]
with \(P_{B_k}(U):=\langle U,B_k\rangle_{x_k}B_k\). A central identity is
\[
\theta=\langle \mathrm{grad}\,f(X),B\rangle_X=\langle \nabla f(X),B\rangle,
\]
which permits coordinate derivatives to be computed directly from the Euclidean gradient [2406.02225].

When \(M=\mathbb{R}^n\), these constructions reduce to classical Euclidean block coordinate descent: tangent spaces become the ambient space, \(\Exp_x(v)=x+v\), and the update becomes \(x^+=x-\alpha P_B\nabla f(x)\) [1912.10627] [2406.02225].

## 2. Update rules and algorithmic families

The product-manifold literature distinguishes three principal block updates. **Block Exact Minimization (BEM)** updates block \(i\) by solving
\[
x_i^{t+1}\in \arg\min_{u\in M_i}F(x_1^{t+1},\dots,x_{i-1}^{t+1},u,x_{i+1}^t,\dots,x_b^t).
\]
**Block Majorization–Minimization (BMM)** replaces the block objective by a block-\(i\) majorizer \(G_i\) and minimizes that surrogate. **Block Riemannian Gradient Descent (BRGD)** uses
\[
x_i^{t+1}\gets \mathrm{Retr}_{x_i^t}(s_i^t),\qquad
s_i^t=-\lambda_i^t\,\widetilde{\nabla}_iF(x_1^{t+1},\dots,x_{i-1}^{t+1},x_i^t,x_{i+1}^t,\dots,x_b^t).
\]
A blended scheme applies BRGD on some blocks and exact minimization on others [2305.14744].

The constrained block-Riemannian optimization literature places these methods in a block majorization-minimization template. At iteration \(k\), a block \(b_k\) is chosen, a surrogate \(Q_{b_k}^k\) is built so that \(Q_{b_k}^k(z)\ge f(x_1^k,\dots,z,\dots,x_B^k)\) and \(Q_{b_k}^k(x_{b_k}^k)=f(x^k)\), and the update minimizes the surrogate over the feasible block set. A first-order quadratic tangent-space surrogate recovers RBCD as a prox-linear special case [2312.10330].

On matrix manifolds, cyclic, random, and Gauss–Southwell selection rules are natural. A fixed-order cyclic sweep updates blocks \(B_1,\dots,B_{|\mathcal J|}\), a random rule samples a block uniformly, and a Gauss–Southwell rule chooses the block maximizing \(|\langle \nabla f(x),B\rangle|\). The same literature also introduces a first-order approximation variant, RCDlin, which anchors the Euclidean gradient at \(X_k\), performs \(S\) inner block updates with that fixed gradient, and evaluates \(\nabla f\) only once per outer iteration [2406.02225].

| Family | Representative update | Setting |
|---|---|---|
| BEM | exact minimization over one block | product manifolds |
| BMM | minimization of a block majorizer \(G_i\) or \(Q_b^k\) | constrained block-Riemannian optimization |
| BRGD | \(x_i^{t+1}=\mathrm{Retr}_{x_i^t}(-\lambda_i^t\widetilde{\nabla}_iF)\) | smooth embedded product manifolds |
| TSD | \(x_{k+1}=R_{x_k}(-\alpha_kP_{S_k}\nabla f(x_k))\) | tangent-subspace blocks |
| RCDlin | \(X_{k}^{s+1}=R_{X_k^s}(-\eta\,\theta_k^s B_{\ell_k^s})\) with anchored \(\nabla f(X_k)\) | matrix manifolds |

The distinction among these schemes is substantive. In BEM and BMM, the block step may be a small exact or surrogate minimization. In BRGD and TSD, the step is explicitly first order. In matrix-manifold RCD, the blocks are often chosen to exploit closed-form retractions with very low per-update cost [2305.14744] [2312.10330] [1912.10627] [2406.02225].

## 3. Geometry of blocks on matrix and manifold domains

The Stiefel manifold,
\[
\mathrm{St}(d,k)=\{U\in\mathbb{R}^{d\times k}\mid U^\top U=I_k\},
\]
is central in RBCD. Its tangent space is
\[
T_U\mathrm{St}(d,k)=\{\Xi\in\mathbb{R}^{d\times k}\mid U^\top\Xi+\Xi^\top U=0\},
\]
and the Euclidean-metric projection is
\[
\mathrm{Proj}_{T_U\mathrm{St}}(Z)=Z-U\,\mathrm{Sym}(U^\top Z).
\]
Common retractions include the QR retraction and the polar retraction. For matrix-manifold coordinate descent on Stiefel, basis directions of the form \(B_{ij}=H_{ij}X\), with \(H_{ij}=e_ie_j^\top-e_je_i^\top\), yield Givens-rotation updates
\[
X^+=G_{ij}(-\eta\theta_{ij})X,
\]
which modify only rows \(i\) and \(j\) and cost \(O(p)\) per update [2012.05199] [2406.02225].

The Grassmann manifold uses the horizontal tangent representation
\[
T_X\mathrm{Gr}(p,n)=\{\Delta\in\mathbb{R}^{n\times p}\mid X^\top\Delta=0\},
\qquad
\mathrm{proj}_X(G)=(I-XX^\top)G.
\]
The same Givens-style update is valid, and the matrix-manifold literature proves that it commutes with right actions, making it well defined on equivalence classes \([X]=\{XQ:Q\in O(p)\}\) [2406.02225].

TSD constructs blocks directly in tangent spaces. On \(O(n)\), the tangent space is \(T_YO_n=\{YA:A\in\mathrm{Skew}_n\}\), and block subspaces are induced by partitions of an orthonormal skew basis into sets \(B_k\). Under a “disjoint pairs” assumption, the exponential of the block generator decomposes into independent \(2\times2\) Givens blocks. On \(\mathrm{St}(n,p)\), the randomized rule uses one-dimensional subspaces spanning skew and normal components, with sampled directions \(UH_{ij}\) or \(ve_{\ell,p}^\top\), where \(v\in\mathrm{Ker}(U^\top)\) is sampled uniformly from the unit sphere [1912.10627].

The matrix-manifold framework extends the same principle to other geometries. For fixed-rank SPSD matrices represented as \(X=YY^\top\), the tangent space at \(Y\) is Euclidean, the gradient is \(\mathrm{grad}\,f(Y)=2\nabla f(YY^\top)Y\), and single-entry updates
\[
Y^+=Y-\eta\theta_{ij}e_ie_j^\top
\]
have \(O(1)\) update cost. For SPD matrices with the Bures–Wasserstein metric, coordinate steps use basis directions \(E_{ij}X+XE_{ij}\) and alter only two rows and columns. For generalized hyperbolic manifolds, symplectic manifolds, doubly stochastic manifolds, and multinomial manifolds, the literature supplies manifold-specific bases, closed-form coordinate derivatives, and retractions such as Lorentz boosts, symplectic exponentials, coordinate Sinkhorn, and row normalization [2406.02225].

A recurrent structural fact is that the block is defined in the tangent geometry, not merely in ambient coordinates. This is why apparently similar “coordinate descent” constructions can differ materially across manifolds.

## 4. Assumptions, convergence guarantees, and complexity regimes

The TSD analysis assumes a geodesically complete Riemannian manifold, smoothness in the form
\[
\|\Gamma_x^y\nabla f(x)-\nabla f(y)\|_y\le L_f d(x,y),
\]
and either the deterministic \((\gamma,r)\)-gap ensuring condition or the randomized \(C\)-randomized norm condition. Under block smoothness and \((\gamma,r)\)-gap ensuring, the deterministic multi-block method satisfies \(\lim_{t\to\infty}\|\nabla f(x^t)\|=0\), any limit point is stationary, the stationarity rate is \(O(m/\sqrt t)\), and for \(g\)-convex objectives the function-value rate is \(O(m^2/t)\). Under the \(C\)-randomized norm condition, the randomized single-block method yields expected stationarity
\[
\min_{s\le t}\mathbb E\|\nabla f(x^{s-1})\|
\le \sqrt{\frac{2L_f(f(x^0)-f^*)}{C^2t}},
\]
almost sure convergence of \(\|\nabla f(x^t)\|\) to zero, and an \(O(1/t)\) expected function-value bound under \(g\)-convexity [1912.10627].

The product-manifold BRGD analysis assumes compact smooth submanifolds, block-\(i\) Lipschitz smoothness of the Euclidean partial gradients, and smooth retractions with blockwise descent lemmas. With cyclic BRGD and stepsizes \(\lambda_i^t=1/L_i\), the theory gives
\[
\min_{t=0,\dots,T}\|\widetilde{\nabla}F(x^t)\|
\le \sqrt{b}\,C_b\,\sqrt{\frac{F(x^0)-F(x^{T+1})}{T+1}},
\]
hence iteration complexity \(\mathcal O(C_b^2(F(x^0)-F^*)/\varepsilon^2)\) to reach \(\varepsilon\)-stationarity. The blended BEM+RGD scheme admits an analogous bound with \(b-1\) gradient-updated blocks [2305.14744].

The RBMM analysis addresses smooth nonconvex objectives over \(X=X_1\times\cdots\times X_B\subseteq M_1\times\cdots\times M_B\), with geodesic smoothness, compact sublevel sets, a uniform injectivity radius bound, and either \(g\)-smooth, Riemannian proximal, or Euclidean proximal surrogates. Every limit point is stationary under the stated assumptions. For proximal surrogates, the worst-case iteration complexity for the paper’s \(\varepsilon\)-stationarity measure is \(\widetilde O(\varepsilon^{-2})\); for \(g\)-smooth surrogates, the general complexity is \(\widetilde O(\varepsilon^{-4})\), improved to \(\widetilde O(\varepsilon^{-2})\) under an upper quadratic surrogate-gap bound. On products of Euclidean and Stiefel manifolds, the assumptions become completely Euclidean while the analysis remains Riemannian [2312.10330].

The matrix-manifold RCD and RCDlin analyses assume a compact neighborhood around a critical point, bounded basis and projection constants, and retraction \(L\)-smoothness. Randomized RCD with one inner step and \(\eta=\Theta(1/L)\) satisfies
\[
\min_{0\le k\le K-1}\|\mathrm{grad}\,f(X_k)\|_{X_k}^2\le O\!\left(\frac{|\mathcal J|\,L}{K}\right),
\]
whereas cyclic RCD yields
\[
\min_{0\le k\le K-1}\|\mathrm{grad}\,f(X_k)\|_{X_k}^2\le O\!\left(\frac{|\mathcal J|^2\,L}{K}\right).
\]
For RCDlin, the randomized rate becomes \(O(|\mathcal J|L/(KS))\) under a positive-correlation condition between anchored and true block derivatives [2406.02225].

A distinct convergence template appears in non-orthogonal joint approximate diagonalization on \(\mathrm{St}(m,n,\mathbb C)\times SL_m(\mathbb C)\). There, gradient-based BCD-G algorithms choose the block according to the Riemannian gradient magnitude, use Wolfe line search on the Stiefel block, and use structured elementary transformations on the special linear block. Under bounded iterates, the whole sequence converges and every limit point is stationary, with the proof based on sufficient descent and the Łojasiewicz gradient inequality [2009.13377].

These results are not interchangeable. Their rates, stationarity measures, and assumptions depend on whether blocks are product components, tangent subspaces, or manifold-specific coordinates.

## 5. The Projection Robust Wasserstein instance

A prominent specialized use of the acronym is the RBCD algorithm for computing the Projection Robust Wasserstein (PRW) distance. For empirical measures \(\mu_n=\sum_{i=1}^n r_i\delta_{x_i}\) and \(\nu_n=\sum_{j=1}^n c_j\delta_{y_j}\), the squared PRW distance is
\[
P_k^2(\mu_n,\nu_n)
=
\max_{U\in \mathrm{St}(d,k)}
\min_{\pi\in\Pi(\mu_n,\nu_n)}
\sum_{i,j=1}^n \pi_{ij}\,\|U^\top x_i-U^\top y_j\|^2.
\]
The entropy-regularized version introduces
\[
H(\pi):=-\sum_{i,j}(\pi_{ij}\log\pi_{ij}-\pi_{ij}),
\]
and maximizes over \(U\in\mathrm{St}(d,k)\) the inner regularized OT value. The key reformulation dualizes the inner RegOT and obtains a smooth minimization problem over two Euclidean blocks \(u,v\in\mathbb R^n\) and one Stiefel block \(U\):
\[
\min_{U\in\mathrm{St}(d,k),\,u,v\in\mathbb R^n}
g(u,v,U)
=
\log\|\zeta(u,v,U)\|_1-r^\top u-c^\top v,
\]
where
\[
[\zeta(u,v,U)]_{ij}
=
\exp\left(-\frac{1}{\eta}\|U^\top(x_i-y_j)\|^2+u_i+v_j\right),
\qquad
\pi=\zeta/\|\zeta\|_1.
\]
This removes the need to solve an entropy-regularized optimal transport problem in each iteration [2012.05199].

The resulting RBCD iteration has three blocks. The \(u\)- and \(v\)-updates are closed-form marginal scaling steps,
\[
u^{t+1}=u^t+\log\big(r\,./\,\varphi(\zeta(u^t,v^t,U^t))\big),
\]
\[
v^{t+1}=v^t+\log\big(c\,./\,\kappa(\zeta(u^{t+1},v^t,U^t))\big),
\]
which are exactly the row and column scaling steps of Sinkhorn applied analytically to \(\zeta\). The \(U\)-block performs one Riemannian gradient step with retraction:
\[
\xi^{t+1}
=
\mathrm{Proj}_{T_{U^t}\mathrm{St}}
\left(
-\frac{2}{\eta}V_{\pi^t}U^t
\right),
\qquad
U^{t+1}=\mathrm{Retr}_{U^t}(-\tau\,\xi^{t+1}),
\]
with
\[
V_\pi=\sum_{i,j}\pi_{ij}(x_i-y_j)(x_i-y_j)^\top.
\]
Feasibility of the transport matrix is enforced at the end by a rounding procedure that corrects marginals and yields \(\hat\pi\in\Pi(\mu_n,\nu_n)\).

The theoretical gain is substantial. The arithmetic complexity to obtain an \((\epsilon_1,\epsilon_2)\)-stationary point is
\[
O\!\left((n^2dk+dk^2+k^3)\log n\left(\frac{1}{\epsilon_2^3}+\frac{1}{\epsilon_1^2\epsilon_2}\right)\right),
\]
and in the common case \(\epsilon_1=\epsilon_2=\epsilon\) and \(k=O(1)\), this becomes \(O((n^2dk+dk^2+k^3)\log n\cdot \epsilon^{-3})\). The comparison paper reports \(O(\epsilon^{-12})\) complexity for RGAS, which solves a Sinkhorn subproblem to high accuracy at every iteration. Per iteration, RBCD costs \(O(n^2dk+dk^2+k^3)\), with \(u\)- and \(v\)-updates of \(O(n)\), \(V_\pi U\) computed in \(O(n^2dk)\) without forming \(V_\pi\), and QR or polar retraction in \(O(dk^2+k^3)\).

The experiments cover synthetic fragmented hypercube and Gaussian models and real data including movie scripts, Shakespeare plays, and MNIST. On the fragmented hypercube with \(k=k^*=2\), the reported runtimes include \(0.18\) s, \(1.83\) s, and \(8.51\) s for RBCD at \(n=100,500,1000\), versus \(0.40\) s, \(2.90\) s, and \(10.25\) s for RGAS. For large scale \(d=250\), \(n=2500\), the reported runtime is \(354.91\) s for RBCD versus \(427.24\) s for RGAS. The paper also reports similar PRW values across methods, robustness to white noise in the Gaussian setting, and an adaptive variant, RABCD, with empirically faster convergence and the same order complexity up to a \(1/\alpha\) factor [2012.05199].

## 6. Applications, relations, and limitations

RBCD appears in a broad range of nonconvex manifold problems. The product-manifold framework covers low-dimensional structure problems such as maximal coding rate reduction, neural collapse, generalized PCA, and alternating projection; combinatorial structure problems such as homomorphic sensing, regression without correspondences, real phase retrieval, and robust point matching; geometric-vision problems such as essential matrix estimation and absolute pose estimation; and outlier-robust estimation via iteratively-reweighted least squares. The same theory recovers previously known results for optimal transport, matrix factorization, and Burer–Monteiro factorization, and yields explicit corollaries for GPCA, essential matrix estimation, PRWD, IRLS, and block-diagonal SDP factorizations [2305.14744].

The matrix-manifold coordinate literature emphasizes computational efficiency in applications such as orthogonal Procrustes on Stiefel, PCA on Grassmann, orthogonal deep network distillation, nearest matrix on the symplectic manifold, Lorentz embeddings, and weighted least squares on SPD or SPSD manifolds. Its practical message is that row-pair rotations, single-entry factor updates, and \(2\times2\) coordinate Sinkhorn updates can be markedly cheaper than full Riemannian gradient or trust-region steps, especially when retraction costs dominate [2406.02225].

In non-orthogonal joint approximate diagonalization, RBCD operates on the product of the complex Stiefel manifold and the special linear group. One block uses Riemannian line-search descent on the Stiefel factor, and the other uses plane Givens, triangular, or diagonal transformations on the special linear factor, producing BCD-GLU, BCD-GQU, and BCD-GU. The square-case Jacobi analogues on \(SL_m(\mathbb C)\) are Jacobi-GLU and Jacobi-GQU. Under bounded iterates, the whole sequence converges and the gradient vanishes [2009.13377].

Several limitations recur across the literature. First, nonconvexity remains intrinsic: the general guarantees are to stationary points or \(\varepsilon\)-stationary points, not to global minimizers. Second, block or subspace selection matters. TSD gives counterexamples showing that poorly chosen subspaces can prevent convergence even if they span the tangent space at each step; the deterministic gap ensuring and randomized \(C\)-randomized norm conditions were introduced precisely to exclude such failures [1912.10627]. Third, low per-iteration cost does not remove all large-scale bottlenecks. In PRW, the \(O(n^2)\) dependence persists even though the method avoids solving a full Sinkhorn problem at every iteration, and the performance remains sensitive to \(\eta\) and \(\tau\) [2012.05199]. Fourth, different RBCD papers use different block notions—product blocks, tangent subspaces, or manifold-specific coordinates—so results from one setting should not be transferred mechanically to another. A plausible implication is that “RBCD” is best understood as a methodological class rather than a single canonical algorithm.

Across these formulations, the common design principle is stable: restrict the search to a geometrically meaningful block, compute a blockwise Riemannian descent direction, and return to the manifold by \(\Exp\) or a retraction. The differences lie in how the blocks are defined, which assumptions are used to prove descent, and how aggressively the algorithm exploits manifold structure for closed-form or low-cost updates.

Source: https://www.emergentmind.com/topics/riemannian-block-coordinate-descent-rbcd