Local convergence of simultaneous min-max algorithms to differential equilibrium on Riemannian manifold

## Abstract

Abstract We study min-max algorithms to solve zero-sum differential games on
Riemannian manifold.
Based on the notions of
differential Stackelberg equilibrium
and differential Nash equilibrium on Riemannian manifold,
we analyze the local convergence of
two representative deterministic simultaneous algorithms τ 𝜏 \tau italic_τ -GDA and τ 𝜏 \tau italic_τ -SGA
to such equilibrium.
Sufficient conditions are obtained to establish their linear convergence rates
by Ostrowski theorem on manifold and spectral analysis.
The τ 𝜏 \tau italic_τ -SGA algorithm is extended from
the symplectic gradient-adjustment method in Euclidean space
to avoid strong rotational dynamics in τ 𝜏 \tau italic_τ -GDA.
In some cases, we obtain a faster convergence rate of τ 𝜏 \tau italic_τ -SGA
through an asymptotic analysis which is valid when
the learning rate ratio τ 𝜏 \tau italic_τ is big.
We show numerically how the insights obtained from the
convergence analysis may improve
the training of orthogonal Wasserstein GANs using
stochastic τ 𝜏 \tau italic_τ -GDA and τ 𝜏 \tau italic_τ -SGA on simple benchmarks.

## 1 Introduction

Riemannian min-max problem has attracted a lot of research attention in recent years,
with various machine learning applications including robust PCA (Jordan et al., ),
robust neural network training (Huang & Gao, )
and generative adversarial network (GAN) (Han et al., ).
This problem is formalized as a two-player zero-sum game,
where the variables of each player
are constrained on the Riemannian manifold $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$,

$$ $\min_{x\in\mathcal{M}_{1}}\max_{y\in\mathcal{M}_{2}}f(x,y).$ $$

When $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean,
gradient-based methods such as
gradient-descent-ascent (GDA) (Daskalakis & Panageas, ; Lin et al., ),
extra-gradient (Gidel et al., ; Mahdavinia et al., ),
optimistic mirror descent (Mertikopoulos et al., ),
Hamiltonian-gradient descent (Loizou et al., )
and symplectic gradient-adjustment (SGA) (Balduzzi et al., )
are often considered to solve this problem.
When $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Riemannian,
how to extend existing algorithms from Euclidean space
to Riemannian manifold
and how to analyze their convergence become
an interesting topic
in recent years.
In Table 1, we summarize
some related works on the Riemannian min-max problem.
In particular, we observe that under suitable assumptions
of the game, one
can obtain the global convergence of a suitable algorithm
towards Nash equilibrium.
One major issue is that these assumptions
typically do not hold in applications such as GAN (Razaviyayn et al., ).
However, when
a min-max problem is non-convex and non-concave,
even in the Euclidean case, the global convergence of
existing algorithms can be very complicated (Hsieh et al., )
and therefore a suitable notion of solution set is needed (Jin et al., ).
In this article, we focus on the notion of differential equilibrium,
and on the local convergence of optimization algorithms
for smooth non-convex non-concave Riemannian min-max problems,
which are relevant to the GAN training.

Section 2
reviews the notion of differential Stackelberg equilibrium (DSE),
which includes the notion of differential Nash equilibrium (DNE)
(Ratliff et al., )
as a special case.
We then analyze the local convergence of
simultaneous min-max algorithms to DSE and DNE,
in which the variables $(x,y)$ are updated simultaneously.
In Section 3.1,
we adopt the classical Ostrowski theorem on manifold
to analyze the local convergence of
deterministic simultaneous
algorithms to a fixed point.
The problem is reduced to analyze the eigenvalues
of a Jacobian matrix in a Euclidean local coordinate,
which is at the heart of analyzing various Euclidean
min-max algorithms (Daskalakis & Panageas, ; Azizian et al., ; Fiez & Ratliff, ; Zhang et al., ; Li et al., ; de Montbrun & Renault, ).
In Section 3.2, sufficient conditions on $\tau$
are given to ensure the local convergence of $\tau$-GDA to DSE or DNE,
where the learning ratio $\tau$ is used to adjust $x$ and $y$ at two different learning rates.

One issue of $\tau$-GDA is its
slow convergence rate
when there are strong rotational dynamics (Li et al., ).
This is a well-known phenomenon near a Nash equilibrium in bilinear games
due to the competitive nature between two players.
To improve the convergence,
first-order methods such as extra-gradient and
optimistic mirror descent are often used to correct the gradient direction of $\tau$-GDA.
Recently, these methods have been extended to Riemannian manifold
(Zhang et al., ; Hu et al., ; Wang et al., ).
However, their computations rely on exponential map and parallel transport
which can be costly (Absil et al., ).
In Section 3.3, we develop an algorithm $\tau$-SGA
to address the same rotational problem based on auto-differentiation.
It naturally extends
the SGA algorithm (Gemp & Mahadevan, ; Balduzzi et al., ; Letcher et al., )
to Riemannian manifold using a learning rate ratio $\tau$.
We provide an asymptotic convergence analysis of deterministic $\tau$-SGA to DSE which is valid when $\tau$ is large.

In Section 4,
we apply $\tau$-GDA and $\tau$-SGA to train orthogonal Wasserstein GANs
(Müller et al., ).
The underlying min-max problem is Riemannian on $\mathcal{M}_{2}$
since we shall impose Stiefel manifold constraints on $y$
to construct Lipschitz-continuous discriminators.
This allows one to compute an approximate Wasserstein distance which can generalize
in high dimension with a polynomial number of training samples (Arora et al., ; Biau et al., ).
Section 5 concludes with some open questions.
In summary, our main contributions are:

$\bullet$
Based on the notions of DSE and DNE on Riemannian manifold,
we derive intrinsic sufficient conditions in terms of the range of $\tau$ and the learning rate of $x$
to obtain a linear convergence rate of $\tau$-GDA to DSE and DNE.

$\bullet$ We develop an algorithm $\tau$-SGA to improve the convergence of deterministic $\tau$-GDA. When $\tau$ is large, we find that a broader range of $\tau$ can be chosen to ensure its convergence to DSE. It indicates that $\tau$-SGA can have a faster local convergence rate compared to $\tau$-GDA.

$\bullet$
We apply the insights from the local convergence analysis to improve the training of orthogonal Wasserstein GANs. We illustrate numerically that an improved convergence of stochastic $\tau$-GDA and $\tau$-SGA may improve the learnt generator on simple benchmarks.(^1^11All the results can be reproduced from a Pytorch software (to be released).)

**Table 1: Related works of Riemannian min-max problem**
| Reference | Class of problem | Solution set / Algorithm |
| --- | --- | --- |
| Zhang et al. () | Geodesically convex (compact) $\mathcal{M}_{1}$,$\mathcal{M}_{2}$ | Nash equilibrium / Extra-gradient |
|  | $f$ is geodesically convex/concave (quasi), semi-continuous (lower/upper) |  |
| Jordan et al. () | +Bounded $\mathcal{M}_{1}$,$\mathcal{M}_{2}$ | Nash equilibrium / Extra-gradient |
|  | $f$ is geodesically convex/concave, smooth |  |
| Huang & Gao () | Euclidean $\mathcal{M}_{2}$ | Stationary point of $\max_{y}f$ / GDA |
|  | $f$ is strongly concave in $y$ |  |
| Han et al. () | Complete $\mathcal{M}_{1}$,$\mathcal{M}_{2}$ | Stationary point of $f$ / Hamiltonian |
|  | $(x,y)\mapsto\|\mbox{grad}_{x}f\|_{x}^{2}+\|\mbox{grad}_{y}f\|_{y}^{2}$ is Polyak-Łojasiewicz |  |
| This work | $f$ is twice<br>continuously differentiable | Differential equilibrium / GDA,SGA |

Notations:
We say that $f\in C^{2}$ if
it is twice continuously differentiable
on the product manifold $\mathcal{M}_{1}\times\mathcal{M}_{2}$.
Let $d_{1}$ denote the dimension of $\mathcal{M}_{1}$,
and $d_{2}$ denote the dimension of $\mathcal{M}_{2}$.
We denote the tangent space of
$\mathcal{M}_{1}$ at $x\in\mathcal{M}_{1}$ by $T_{x}\mathcal{M}_{1}$,
and the tangent space of $\mathcal{M}_{2}$ at $y\in\mathcal{M}_{2}$ by $T_{y}\mathcal{M}_{2}$.
Let $D_{x}$ and $\nabla_{x}$ (resp. $D_{y}$ and $\nabla_{y}$) denote
the differential operator and the Riemannian connection
on $\mathcal{M}_{1}$ (resp. on $\mathcal{M}_{2}$).
Let $\partial$ and $\partial^{2}$ denote the first-order and
second-order partial derivatives of a function on Euclidean space.
For a real symmetric matrix $A$,
$\lambda_{max}(A)$ denotes the maximal eigenvalue of $A$,
$\lambda_{min}(A)$ denotes the minimal eigenvalue of $A$,
and $\lambda_{k}(A)$ denotes its $k$-th smallest eigenvalue.
For a linear transform $B$, $\rho(B)$ denotes its spectral radius
and $\|B\|$ denotes its operator norm.

## 2 Differential equilibrium on Riemannian manifold

This section reviews the definition of DSE,
which includes DNE as a special case.
When $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean,
these solution sets are generic
among local Stackelberg equilibrium (resp. local Nash equilibrium)
of smooth $f$ (Fiez et al., ).
We provide some simple examples to illustrate their difference.

### 2.1 Differential Stackelberg equilibrium (DSE)

When $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean,
the notion of DSE is defined based
on the first and second order derivatives of $f$
at this equilibrium point (Fiez et al., ).
The next definition of DSE is a natural extension
of this concept to Riemannian manifold.

###### Definition 2.1 .

We say that $(x^{\ast},y^{\ast})$ is a DSE of $f\in C^{2}$ if

$$ $\displaystyle\mbox{grad}_{x}f(x^{\ast},y^{\ast})=0,\quad\mbox{grad}_{y}f(x^{ \ast},y^{\ast})=0,$ (1) $\displaystyle-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})\quad\mbox{p.d},$ (2) $\displaystyle[\mbox{Hess}_{x}f-\mbox{grad}^{2}_{yx}f\cdot(\mbox{Hess}_{y}f)^{- 1}\cdot\mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast})\quad\mbox{p.d}.$ (3) $$

In this definition, we rely on the following intrinsic concepts
in the literature of Riemannian optimization (see e.g. Han et al. ()):

$\bullet$ Riemannian gradient: $\mbox{grad}_{x}f(x,y)\in T_{x}\mathcal{M}_{1},\mbox{grad}_{x}f(x,y)\in T_{y}
\mathcal{M}_{2}$,
$\bullet$ Riemannian Hessian: $\mbox{Hess}_{x}f(x,y):T_{x}\mathcal{M}_{1}\to T_{x}\mathcal{M}_{1},\mbox{Hess}
_{y}f(x,y):T_{y}\mathcal{M}_{2}\to T_{y}\mathcal{M}_{2}$,
$\bullet$ Riemannian cross-gradient: $\mbox{grad}^{2}_{xy}f(x,y):T_{x}\mathcal{M}_{1}\to T_{y}\mathcal{M}_{2},\mbox{
grad}^{2}_{yx}f(x,y):T_{y}\mathcal{M}_{2}\to T_{x}\mathcal{M}_{1}$.

The condition (1) means that $(x^{\ast},y^{\ast})$ is a critical point of $f$. Note that the eigenvalues of $\mbox{Hess}_{x}f(x^{\ast},y^{\ast})$ and $\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ depend on the Riemannian metric on $\mathcal{M}_{1}\times\mathcal{M}_{2}$.
However, the notion of DSE does not depend on the choice of Riemannian metric, due to the following known fact.

Fact: The notion of DSE in Definition 2.1
can be equivalently defined in a local coordinate chart which
does not depend on the choice of Riemannian metric.

To make this point clear, we take a local coordinate chart
$(O_{1}\times O_{2},\varphi_{1}\times\varphi_{2})$
around $(x^{\ast},y^{\ast})\in\mathcal{M}_{1}\times\mathcal{M}_{2}$.
We rewrite the function $f(x,y)$ on $(x,y)\in O_{1}\times O_{2}$ using
this chart by

$$ $\bar{f}(u_{1},u_{2})=f(\varphi^{-1}_{1}(u_{1}),\varphi^{-1}_{2}(u_{2})).$ $$

In Appendix A, we verify
the equivalence
between (1)-(3) and the following conditions
(4)-(6):

$$ $\displaystyle\partial_{u_{1}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})=0,\quad \partial_{u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})=0,$ (4) $\displaystyle-\partial_{u_{2}u_{2}}^{2}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})\quad \mbox{p.d (abbrev. of positive definite)},$ (5) $\displaystyle\left[\partial_{u_{1}u_{1}}^{2}\bar{f}-\partial_{u_{2}u_{1}}^{2} \bar{f}\cdot(\partial_{u_{2}u_{2}}^{2}\bar{f})^{-1}\cdot\partial_{u_{1}u_{2}}^ {2}\bar{f}\right](u_{1}^{\ast},u_{2}^{\ast})\quad\mbox{p.d}.$ (6) $$

We see that the conditions (4)-(6) do not depend
on the choice of Riemannian metric,
therefore Definition 2.1 still holds
if the Riemannian metric is changed on $\mathcal{M}_{1}\times\mathcal{M}_{2}$.

It is known that a DSE is a local minimax point (Jin et al., ).
These conditions imply that
$(u_{1}^{\ast},u_{2}^{\ast})=(\varphi_{1}(x^{\ast}),\varphi_{2}(y^{\ast}))$ is a local minimax point of $\bar{f}$. Furthermore, from
the implicit function theorem on Riemannian manifold in Appendix B,
$(x^{\ast},y^{\ast})$ is a local minimax point of $f$ in the following sense:
there exists an open subset $U_{1}\times U_{2}$ of $\mathcal{M}_{1}\times\mathcal{M}_{2}$,
which includes $(x^{\ast},y^{\ast})$ and on which there is a unique function $h:U_{1}\to U_{2}$,
such that the following holds

$$ $f(x^{\ast},y)\leq f(x^{\ast},y^{\ast})\leq\max_{y^{\prime}\in U_{2}}f(x,y^{ \prime})=f(x,h(x)),\quad\forall(x,y)\in U_{1}\times U_{2}.$ $$

The set $\{(x,h(x))|x\in U_{1}\}$ is sometimes called the ridge near the DSE (Wang et al., ).

###### Definition 2.1 .

### 2.2 Differential Nash equilibrium (DNE) and examples

In game theory, one is often interested in finding
a Nash equilibrium since it maintains a symmetry between the role of the players.
When $\mathcal{M}_{1}$
and $\mathcal{M}_{2}$ are Euclidean, it is also called
“strongly local min-max point”
(Daskalakis & Panageas, , Definition 1.6).

The notion of DNE was introduced in Ratliff et al. ()[Definition 3]
through a local coordinate chart. This is equivalent
to the following intrinsic definition:

###### Definition 2.2 .

We say that $(x^{\ast},y^{\ast})$ is a DNE of $f\in C^{2}$ if

$$ $\displaystyle\mbox{grad}_{x}f(x^{\ast},y^{\ast})=0,\quad\mbox{grad}_{y}f(x^{ \ast},y^{\ast})=0,$ (7) $\displaystyle-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})\quad\mbox{p.d},\quad\mbox{ Hess}_{x}f(x^{\ast},y^{\ast})\quad\mbox{p.d}.$ (8) $$

From the definition, it is clear that a DNE is a DSE.
We remark that this concept is defined locally,
and therefore it is different to global Nash-type equilibria on manifold
(Kristály, ).

###### Definition 2.2 .

#### 2.2.1 Example 1: DSE

Consider $f(x,y)=\langle y,Ax-b\rangle$,
with $x\in\mathcal{M}_{1}=\mathbb{R}^{d_{1}}$
and $y\in\mathcal{M}_{2}=S^{d_{2}}$.
The manifold $S^{d_{2}}$ is the unit sphere
embedded in $\mathbb{R}^{d_{2}+1}$, endowed with
the Euclidean metric $\langle\cdot,\cdot\rangle$ on $\mathbb{R}^{d_{2}+1}$.
Let $A^{+}$ denote the pseudo-inverse of $A\in\mathbb{R}^{(d_{2}+1)\times d_{1}}$.
We next present a senario where DSE exists.

###### Proposition 2.1 .

Assume $b\not\in\mbox{Range}(A)$, $\mbox{Ker}(A)=\{0\}$.
Let $x^{\ast}=A^{+}b$, $y^{\ast}=\frac{Ax^{\ast}-b}{\|Ax^{\ast}-b\|}$,
then $(x^{\ast},y^{\ast})$ is a DSE of the $f$ in Example 1.

The proof is given in Appendix C.
Since $\mathcal{M}_{1}$ is Euclidean, it is clear that $\partial^{2}_{xx}f(x^{\ast},y^{\ast})=0$.
This implies that the $(x^{\ast},y^{\ast})$ is not DNE.
But each eigenvalue of $\mathbf{A}=-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ equals to $\|Ax^{\ast}-b\|>0$,
since from the proof $\langle\mathbf{A}[\eta^{\ast}],\eta^{\ast}\rangle=\|Ax^{\ast}-b\|\|\eta^{\ast}
\|^{2}$ for any $\eta^{\ast}\in T_{y^{\ast}}\mathcal{M}_{2}$.

###### Proposition 2.1 .

#### 2.2.2 Example 2: DSE

Consider $f(x,y)=\langle y,Ax-b\rangle-\frac{\kappa}{2}\|Ax\|^{2}$,
with $x\in\mathcal{M}_{1}=\mathbb{R}^{d_{1}}$
and $y\in\mathcal{M}_{2}=S^{d_{2}}$.
Compared to the $f$ in Example 1,
we add a quadratic function of $Ax$ with a curvature parameter $\kappa>0$.
We next show that if $\kappa$ is close to zero, we can still find a DSE near the DSE in
Proposition 2.1.

###### Proposition 2.2 .

Assume $b\not\in\mbox{Range}(A)$, $\mbox{Ker}(A)=\{0\}$.
There exists $\kappa_{0}>0$ such that for any $0<\kappa<\kappa_{0}$,
there is a number $c$ close to $1$,
$x^{\ast}=cA^{+}b$ and $y^{\ast}=\frac{Ax^{\ast}-b}{\|Ax^{\ast}-b\|}$,
so that $(x^{\ast},y^{\ast})$ is a DSE of the $f$ in Example 2.

The proof is given in Appendix D.
From the proof, we have $\mathbf{C}=\partial^{2}_{xx}f(x,y)=-\kappa A^{\intercal}A$.
The spectral radius of $\mathbf{C}$ equals to its operator norm $\|\mathbf{C}\|=\kappa\|A^{\intercal}A\|>0$.
As in Example 1, all the eigenvalues of $\mathbf{A}=-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$
equal to $\|Ax^{\ast}-b\|=\|cAA^{+}b-b\|$.
These quantities will be used in Section 3
to analyze the local convergence of min-max algorithms.

###### Proposition 2.2 .

#### 2.2.3 Example 3: DNE

Consider $f(x,y)=\frac{1}{2}\|Ax+y-b\|^{2}$
with $x\in\mathcal{M}_{1}=\mathbb{R}^{d_{1}}$
and $y\in\mathcal{M}_{2}=S^{d_{2}}$.
$\mathcal{M}_{2}$ is the same embedded sub-manifold of $\mathbb{R}^{d_{2}+1}$ as above.
We next provide a sufficient condition on the existence of DNE. The proof is given in Appendix E.

###### Proposition 2.3 .

Assume $b\not\in\mbox{Range}(A)$, $\mbox{Ker}(A)=\{0\}$.
Let $x^{\ast}=A^{+}b$, $y^{\ast}=\frac{Ax^{\ast}-b}{\|Ax^{\ast}-b\|}$,
then $(x^{\ast},y^{\ast})$ is a DNE of the $f$ in Example 3.

###### Proposition 2.3 .

## 3 Simultaneous Min-max algorithms for differential equilibrium

Simultaneous gradient-based min-max algorithms such as
GDA and SGA are often used
to find local Nash equilibrium
when $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are
Euclidean (Daskalakis & Panageas, ; Letcher et al., ).
Similar to GDA, we extend the SGA algorithm
to Riemannian manifold with two-time scale update (Heusel et al., )
using either deterministic or stochastic gradients.
Section 3.1 reviews a classical result
of fixed point theorem. Based on this theorem,
we then focus on a deterministic analysis of the local convergence of these algorithms
to DSE and DNE.

### 3.1 Local convergence of deterministic simultaneous algorithms

We use the classical Ostrowski theorem
to analyze the local convergence of simultaneous deterministic algorithms on manifold.
Each algorithm is defined by an update rule which does not change over iteration.
This theorem provides a sufficient condition
for the linear convergence rate of an algorithm to a fixed point,
based on the spectral radius of
the update rule’s Jacobian matrix at the fixed point.

A deterministic simultaneous algorithm is defined by
two vector fields $(x,y)\mapsto\xi_{1}(x,y)\in T_{x}\mathcal{M}_{1}$,
$(x,y)\mapsto\xi_{2}(x,y)\in T_{y}\mathcal{M}_{2}$ on $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$,
and a suitable choice of manifold retractions.
Initialized at a point $(x(0),y(0))$ on $\mathcal{M}_{1}\times\mathcal{M}_{2}$,
the algorithm generates a sequence $(x(t+1),y(t+1))=\mathbf{T}(x(t),y(t))$,
through an update rule $\mathbf{T}:\mathcal{M}_{1}\times\mathcal{M}_{2}\to\mathcal{M}_{1}\times
\mathcal{M}_{2}$ of the following form,

$$ $\mathbf{T}(x,y)=(\mathcal{R}_{1,x}(\xi_{1}(x,y)),\mathcal{R}_{2,y}(\xi_{2}(x,y ))),$ $$

where $\mathcal{R}_{1,x}:T_{x}\mathcal{M}_{1}\to\mathcal{M}_{1}$
(resp. $\mathcal{R}_{2,y}:T_{y}\mathcal{M}_{2}\to\mathcal{M}_{2}$)
denotes the restriction of a retraction $\mathcal{R}_{1}$ at $x\in\mathcal{M}_{1}$
(resp. retraction of $\mathcal{R}_{2}$ at $y\in\mathcal{M}_{2}$).
For example, on the
$\mathcal{M}_{1}=\mathbb{R}^{d_{1}}$ and $\mathcal{M}_{2}=S^{d_{2}}$
in Section 2.2,
we will take $\mathcal{R}_{1,x}(\delta)=x+\delta$
for $\delta\in\mathbb{R}^{d_{1}}$,
and $\mathcal{R}_{2,y}(\eta)$ to be the projection of
the vector $y+\eta$ in $\mathbb{R}^{d_{2}+1}$
to the sphere $\mathcal{M}_{2}$.

We say that $(x^{\ast},y^{\ast})\in\mathcal{M}_{1}\times\mathcal{M}_{2}$ is a fixed point of $\mathbf{T}$
if it is a critical point of the vector fields $\xi_{1}$ and $\xi_{2}$.
We next define what it means for
an update rule $\mathbf{T}$ to be locally convergent to
$(x^{\ast},y^{\ast})$.
It implies that it is also a point of attraction of $\mathbf{T}$ (Ortega & Rheinboldt, , Definition 10.1.1).

###### Definition 3.1 (locally convergent with rate $\rho\in(0,1)$ ) .

Let $(x^{\ast},y^{\ast})$ be a fixed point of $\mathbf{T}$.
For any $\epsilon\in(0,1-\rho)$,
there exists a local stable region $S_{\delta}\subset\mathcal{M}_{1}\times\mathcal{M}_{2}$, which is a
geodesically convex open set (and homeomorphic to a ball) containing $(x^{\ast},y^{\ast})$ such that started from $(x(0),y(0))\in S_{\delta}$, we have $(x(t),y(t))\in S_{\delta},\forall t\geq 1$.
Furthermore, let $d(t)$ be the Riemannian distance between $(x(t),y(t))$
and $(x^{\ast},y^{\ast})$, then there exists
a constant $C$ s.t. $d(t)\leq C(\rho+\epsilon)^{t+1}d(0),\forall t\geq 0$.

We next give a sufficient condition of $\mathbf{T}$ to be locally convergent to $(x^{\ast},y^{\ast})$. It is based on the spectral radius of the following linear transformation on $T_{x^{\ast}}\mathcal{M}_{1}\times T_{y^{\ast}}\mathcal{M}_{2}$,

$$ $\mathbf{D}\mathbf{T}^{\ast}=\mathbf{T}^{\prime}(x^{\ast},y^{\ast})=\left( \begin{array}[]{ll}I+\nabla_{x}\xi_{1}(x^{\ast},y^{\ast})&D_{y}\xi_{1}(x^{\ast },y^{\ast})\\ D_{x}\xi_{2}(x^{\ast},y^{\ast})&I+\nabla_{y}\xi_{2}(x^{\ast},y^{\ast})\end{ array}\right).$ (9) $$

###### Theorem 3.1 (Ostrowski Theorem on manifold) .

Let $(x^{\ast},y^{\ast})$ be a fixed point of $\mathbf{T}$.
Assume that $\xi_{1}$ and $\xi_{2}$ are continuous on $\mathcal{M}_{1}\times\mathcal{M}_{2}$,
and they are differentiable at $(x^{\ast},y^{\ast})$
such that $\rho(\mathbf{D}\mathbf{T}^{\ast})<1$,
then $\mathbf{T}$ is locally convergent to $(x^{\ast},y^{\ast})$
with rate $\rho(\mathbf{D}\mathbf{T}^{\ast})$.

This result is proved in Ortega & Rheinboldt (, Section 10.1.3) when
$\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean.
The proof idea can be readily extended to the manifold case.
In the statement of Theorem 3.1, we add an assumption
on the continuity to $\xi_{1}$ and $\xi_{2}$
to construct a non-empty local stable region $S_{\delta}$.
To make the article self-contained, we provide a proof of
this theorem in Appendix F.

Theorem 3.1
has a local nature as
the convergence rate $\rho(\mathbf{D}\mathbf{T}^{\ast})$
does not depend on any global manifold property.
One can also use
the operator norm condition $\|\mathbf{D}\mathbf{T}^{\ast}\|<1$
on the tangent map $\mathbf{D}\mathbf{T}^{\ast}$ of $\mathbf{T}$ at $(x^{\ast},y^{\ast})$
in Boumal (, Theorem 4.19)
to obtain a similar local convergence result.
In this article, we follow Theorem 3.1
as $\rho(\mathbf{D}\mathbf{T}^{\ast})<1$ is a weaker assumption.

###### Definition 3.1 (locally convergent with rate ρ ∈ ( 0 , 1 ) 𝜌 0 1 \rho\in(0,1) italic_ρ ∈ ( 0 , 1 ) ) .

###### Theorem 3.1 (Ostrowski Theorem on manifold) .

### 3.2 Simultaneous gradient-descent-ascent algorithm ( τ 𝜏 \tau italic_τ -GDA)

The $\tau$-GDA algorithm uses
the Riemannian gradients $\mbox{grad}_{x}f(x,y)$ and $\mbox{grad}_{y}f(x,y)$
to update $x$ and $y$ simultaneously.
The local convergence of deterministic $\tau$-GDA to DSE and DNE
is well-studied in Euclidean space (Daskalakis & Panageas, ; Jin et al., ; Fiez & Ratliff, ; Li et al., ).
This section extends these results to Riemannian manifold.
We obtain a sharp lower-bound of $\tau$ for $\tau$-GDA to be locally convergent to DSE.
It is based on a refinement of the spectral analysis
in Euclidean space (Li et al., ).

##### τ 𝜏 \tau italic_τ -GDA algorithm

In deterministic setting,
the update rule $\mathbf{T}$ of $\tau$-GDA is determined by

$$ $\displaystyle\xi_{1}(x,y)=-\gamma\mbox{grad}_{x}f(x,y),\quad\xi_{2}(x,y)=\tau \gamma\mbox{grad}_{y}f(x,y).$ (10) $$

Note that $\gamma>0$ and we use a ratio $\tau>0$ to adjust
the learning rate (step size) of the Riemannian gradients.
The deterministic $\tau$-GDA can be readily extended to stochastic $\tau$-GDA
by using an unbiased estimation of the Riemannian gradients (Jordan et al., ; Huang & Gao, ).

##### Local convergence of deterministic τ 𝜏 \tau italic_τ -GDA

Based on Theorem 3.1,
we are ready to study the local convergence
of $\tau$-GDA (defined by (10))
to DSE and DNE.
We rewrite the Riemannian Hessian and cross-gradient
in the following linear transform on $T_{x^{\ast}}\mathcal{M}_{1}\times T_{y^{\ast}}\mathcal{M}_{2}$ using $(\mathbf{A},\mathbf{B},\mathbf{C})$, i.e.

$$ $\mathbf{M}_{g}=\left(\begin{array}[]{ll}-\mathbf{C}&-\mathbf{B}\\ \tau\mathbf{B}^{\intercal}&-\tau\mathbf{A}\end{array}\right)=\left(\begin{ array}[]{ll}-\mbox{Hess}_{x}f(x^{\ast},y^{\ast})&-\mbox{grad}^{2}_{yx}f(x^{ \ast},y^{\ast})\\ \tau\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})&\tau\mbox{Hess}_{y}f(x^{\ast},y^{ \ast})\end{array}\right).$ (11) $$

For a linear transform $\mathbf{M}$ with only negative eigenvalues, we write
$\gamma^{\bullet}(\mathbf{M})=-2\max_{k}\frac{\mbox{Re}(\lambda_{k}(\mathbf{M})
)}{|\lambda_{k}(\mathbf{M})|}$.
Let $L_{g}=\max(\|\mathbf{A}\|,\|\mathbf{B}\|,\|\mathbf{C}\|)$
and $\mu_{g}=\min(L_{g},\lambda_{min}(\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B
}^{\intercal}))$.

###### Theorem 3.2 .

Assume $(x^{\ast},y^{\ast})$ is a DSE of $f\in C^{2}$.
If $\tau>\frac{\|\mathbf{C}\|}{\lambda_{min}(\mathbf{A})}$ and $\gamma\in(0,\gamma^{\bullet}(\mathbf{M}_{g}))$,
$\tau$-GDA
is locally convergent to $(x^{\ast},y^{\ast})$
with rate $\rho(I+\gamma\mathbf{M}_{g})$.
Furthermore, if $\tau\geq\frac{2L_{g}}{\lambda_{min}(\mathbf{A})}$
and $\gamma=\frac{1}{4\tau L_{g}}$, the rate
is at most $1-\frac{\mu_{g}}{16\tau L_{g}}$.

The proof is given in Appendix G.
This result is an extension of the Euclidean space result
in Li et al. (, Theorem 4.2)
(without assuming $\mathbf{M}_{g}$ being diagonalizable).
When $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean,
the range $\{\tau\in\mathbb{R}_{+}|\tau>\|\partial_{xx}^{2}f(x^{\ast},y^{\ast})\|/\lambda
_{min}(-\partial_{yy}^{2}f(x^{\ast},y^{\ast}))\}$ in Theorem 3.2 is sharp as one can construct
a counter-example
(Li et al., , Theorem 4.1)
to show that if $\tau$ is outside this range,
the spectral radius of
the Jacobian matrix (9)
is strictly larger than one for any $\gamma>0$
(see also a discussion in Li et al. (, Remark 2)).

Theorem 3.2 can be readily applied
to analyze the local convergence of $\tau$-GDA to DNE.
However, we can obtain a broader range of $\tau$.
In the Euclidean case,
the local convergence to DNE is analyzed
in Daskalakis & Panageas (, Theorem 1.10) when $\tau=1$,
and generalized to $\tau>0$ (Jin et al., , Proposition 26).

###### Theorem 3.3 .

Assume $(x^{\ast},y^{\ast})$ is a DNE of $f\in C^{2}$.
If $\tau>0$ and $\gamma\in(0,\gamma^{\bullet}(\mathbf{M}_{g}))$,
$\tau$-GDA is locally convergent to
$(x^{\ast},y^{\ast})$
with rate $\rho(I+\gamma\mathbf{M}_{g})$.

The proof is given in Appendix H.
Different to the proof of Theorem 3.2,
we analyze the spectral radius of
the Jacobian matrix (9)
based on the Ky Fan inequality (Zhan, , Lemma 4.20)
as in Daskalakis & Panageas (, Lemma 2.7).

###### Theorem 3.2 .

###### Theorem 3.3 .

### 3.3 Symplectic gradient-adjustment method ( τ 𝜏 \tau italic_τ -SGA)

In Euclidean space, SGA modifies
the update rule of GDA
to avoid strong rotational dynamics near
a fixed point (Gemp & Mahadevan, ).
We apply this idea to $\tau$-GDA
and extend it to Riemannian manifold.
Its local convergence to DSE is analyzed when $\tau$ is large.

##### τ 𝜏 \tau italic_τ -SGA algorithm

SGA adjusts
a vector field $\xi$ using
an orthogonal vector field which is
constructed from the anti-symmetric part of the
Jacobian matrix of $\xi$.
In $\tau$-GDA, $\xi(x,y)=(-\delta(x,y),\tau\eta(x,y))$
where $\delta(x,y)=\mbox{grad}_{x}f(x,y)$ and
$\eta(x,y)=\mbox{grad}_{y}f(x,y)$.
Based on this idea,
the adjustment of $\tau$-SGA
depends on the Riemannian cross-gradients $\tilde{\mathbf{B}}(x,y)=\mbox{grad}^{2}_{yx}f(x,y)$
and $\tilde{\mathbf{B}}^{\intercal}(x,y)=\mbox{grad}^{2}_{xy}f(x,y)$,
which is summarized in the following update rule (with a hyper-parameter $\mu\in\mathbb{R}$),

$$ $\displaystyle\xi_{1}(x,y)$ $\displaystyle=-\gamma\left(\delta(x,y)+\mu\frac{(\tau+1)\tau}{2}\tilde{\mathbf {B}}(x,y)[\eta(x,y)]\right),$ (12) $\displaystyle\xi_{2}(x,y)$ $\displaystyle=\gamma\left(\tau\eta(x,y)-\mu\frac{\tau+1}{2}\tilde{\mathbf{B}}^ {\intercal}(x,y)[\delta(x,y)]\right).$ (13) $$

In Appendix I, we detail the derivation
of the update rule. The next proposition shows that
the adjusted direction $(-\tau\tilde{\mathbf{B}}[\eta],-\tilde{\mathbf{B}}^{\intercal}[\delta])$
is orthogonal to $\xi=(-\delta,\tau\eta)$.

###### Proposition 3.1 .

For any $(x,y)\in\mathcal{M}_{1}\times\mathcal{M}_{2}$, we have

$$ $\langle\tau\tilde{\mathbf{B}}(x,y)[\eta(x,y)],\delta(x,y)\rangle_{x}+\langle \tilde{\mathbf{B}}^{\intercal}(x,y)[\delta(x,y)],-\tau\eta(x,y)\rangle_{y}=0.$ $$

The proof is given in Appendix I.1.
In the original SGA method, the orthogonality is essential to make
it compatible to potential and Hamiltonian game dynamics. This proposition
implies that such compatibility still makes sense for $\tau$-SGA.
In Appendix K, we discuss how to perform deterministic and
stochastic gradient adjustment
using auto-differentiation when $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$ are Euclidean embedded sub-manifolds
(Absil et al., , Chapter 3.3).

###### Proposition 3.1 .

##### Local convergence of deterministic τ 𝜏 \tau italic_τ -SGA to DSE: asymptotic analysis

We study the local convergence of $\tau$-SGA with deterministic gradients
in an asymptotic regime where $\tau\to\infty$.
In this regime, to make the term $\tilde{\mathbf{B}}[\eta]$
comparable to the term $\delta$ in (12),
we re-parameterize $\mu=\theta\frac{2}{\tau(\tau+1)}\sim\frac{1}{\tau^{2}}$.
As a consequence,
$\mu\frac{\tau+1}{2}\sim\frac{1}{\tau}$
and the update rule in (13) can be
approximated by the $\xi_{2}$ in (10).

The next theorem analyzes the local convergence of Asymptotic $\tau$-SGA,
which is an approximation of $\tau$-SGA,
whose $\xi_{1}$ (resp. $\xi_{2}$) is defined by (12)
(resp. (10)).
In order to apply Theorem 3.1,
it is necessary to verify that the corresponding
$\xi_{1}$ is differentiable at $(x^{\ast},y^{\ast})$.
This is indeed true since $\tilde{\mathbf{B}}(x,y)$
is continuous at $(x^{\ast},y^{\ast})$
and $\eta(x^{\ast},y^{\ast})=0$.
The following linear transform plays
a key role in the asymptotic analysis,

$$ $\mathbf{M}_{s}=\left(\begin{array}[]{ll}-\mathbf{C}&-\mathbf{B}\\ \tau\mathbf{B}^{\intercal}&-\tau\mathbf{A}\end{array}\right)+\theta\left( \begin{array}[]{ll}-\mathbf{B}\mathbf{B}^{\intercal}&\mathbf{B}\mathbf{A}\\ 0&0\end{array}\right).$ (14) $$

Let $L_{s}=\max(\|\mathbf{A}\|,\|\mathbf{B}\|,\|\mathbf{C}+\theta\mathbf{B}\mathbf{
B}^{\intercal}\|)$
and $\mu_{s}=\min(L_{s},\lambda_{min}(\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B
}^{\intercal}))$.

###### Theorem 3.4 .

Assume $(x^{\ast},y^{\ast})$ is a DSE of $f\in C^{2}$.
If $\tau>\min\left(\frac{\|\mathbf{C}\||}{\lambda_{min}(\mathbf{A})},\frac{\|
\mathbf{C}+\theta\mathbf{B}\mathbf{B}^{\intercal}\|}{\lambda_{min}(\mathbf{A})
}\right)$,
$\mu=\theta\frac{2}{\tau(\tau+1)}$ with $0\leq\theta\leq\frac{1}{\lambda_{max}(\mathbf{A})}$,
and $\gamma\in(0,\gamma^{\bullet}(\mathbf{M}_{s}))$,
Asymptotic $\tau$-SGA
is locally convergent to $(x^{\ast},y^{\ast})$
with rate $\rho(I+\gamma\mathbf{M}_{s})$.
Furthermore, if $\tau\geq\frac{2L_{s}}{\lambda_{min}(\mathbf{A})}$
and $\gamma=\frac{1}{4\tau L_{s}}$, the rate
is at most $1-\frac{\mu_{s}}{16\tau L_{s}}$.

The proof is given in Appendix J.
Contrary to Theorem 3.2,
the valid range of $\tau$ depends also on the choice of $\theta$.
Indeed, if the spectral radius of $\mathbf{C}$ is larger than that of $\mathbf{C}+\theta\mathbf{B}\mathbf{B}^{\intercal}$
with a suitable choice of $\theta$,
a broader range of $\tau$ could be used in $\tau$-SGA compared to $\tau$-GDA.
Furthermore,
we obtain a non-trivial improvement on the convergence rate
when $L_{s}=\mu_{s}=\|\mathbf{C}+\theta\mathbf{B}\mathbf{B}^{\intercal}\|<L_{g}=\mu_
{g}=\|\mathbf{C}\|$.
In this case, we could choose a smaller $\tau_{s}=\frac{2L_{s}}{\lambda_{min}(\mathbf{A})}$
and a larger $\gamma_{s}=\frac{1}{4\tau_{s}L_{s}}=\frac{\lambda_{min}(\mathbf{A})}{8L_{s}^{2}}$ in $\tau$-SGA
to achieve an faster rate $1-\frac{1}{16\tau_{s}}$,
compared to the rate $1-\frac{1}{16\tau_{g}}$ using a larger $\tau_{g}=\frac{2L_{g}}{\lambda_{min}(\mathbf{A})}$ and a smaller $\gamma_{g}=\frac{1}{4\tau_{g}L_{g}}$ in $\tau$-GDA.
This rate of Asymptotic $\tau$-SGA is non-trivial
as it significantly improves the rate of the extra-gradient method
in Euclidean space (Li et al., , Theorem 5.4).

We next illustrate the convergence results using a numerical example,
to show that $\tau$-SGA can indeed converge much faster than
$\tau$-GDA to DSE when there are strong rotational forces in its dynamics.
In Figure 1,
we compare
the local convergence rate of $\tau$-GDA and
$\tau$-SGA in Example 2 of Section 2.2,
where $A=[1;1;1]\in\mathbb{R}^{3\times 1}$, $b=[1;1;0.99]\in\mathbb{R}^{3}$
and $\kappa=0.1$.
In this case, we find numerically that $x^{\ast}\approx 0.9975$.
The initial point of each algorithm
is $x=A^{+}b$,$y=(Ax-b)/\|Ax-b\|$.
In Figure 1(a),
we study the range of $\tau$ for $\tau$-GDA convergence.
According to Theorem 3.2
and Proposition 2.2,
a valid range of $\tau$ should be larger than
$\frac{\max_{k}|\lambda_{k}(\mathbf{C})|}{\lambda_{min}(\mathbf{A})}=\kappa\|A^
{\intercal}A\|/\|Ax^{\ast}-b\|\approx 36.18$.
Figure 1(a)
shows that when $\tau=30$,
$\tau$-GDA can slowly diverge. When $\tau=50$, $\tau$-GDA converges
slowly to the value $\max_{y}f(x^{\ast},y)=-0.141$.
This is a convergence dilemma of $\tau$-GDA:
to achieve local convergence, $\tau$ needs to be large,
but the rate can be very slow.
In Figure 1(b) and (c),
we study $\tau$-SGA using the same learning rate $\tau\gamma$ for $y$.
Figure 1(b) shows
that the gap between
Asymptotic $\tau$-SGA
and $\tau$-SGA is not so big even if $\tau=10$.
Therefore Theorem 3.4 is valid for $\tau$-SGA.
When $\theta=0.15$,  we find that
$\frac{\max_{k}|\lambda_{k}(\mathbf{C}+\theta\mathbf{B}\mathbf{B}^{\intercal})|
)|}{\lambda_{min}(\mathbf{A})}\approx 16.7$,
which suggests that $\tau$-SGA can be locally convergent
with a smaller $\tau$ compared to $\tau$-GDA.
Figure 1(c) confirms
this analysis and it also shows that a faster rate can be achieved by $\tau$-SGA.
Surprisingly, $\tau$-SGA quickly converges towards $-0.141$ even if $\tau=10$.

Figure: (a)
Refer to caption: extracted/5893895/plots/iclr_compute_dse_dert_gdatau_example2_tau30.png

The local convergence of $\tau$-SGA
to DNE in Euclidean space is analyzed in
Letcher et al. (, Theorem 10) when $\tau=1$.
However, we could not easily extend this result
to $\tau>0$,
since when $\tau\neq 1$, a DNE is not necessarily
a stable fixed point of the
vector field $\xi=(-\delta,\tau\eta)$,
c.f Letcher et al. (, Definition 4).
It is thus unclear
if there is a non-trivial range of $\mu$
such that for any $\tau>0$, $\tau$-SGA
is locally convergent to DNE.

###### Theorem 3.4 .

## 4 Application to Wasserstein GAN

The local convergence analysis
in Section 3
shows that a larger $\tau$ is sometimes
needed to ensure the local convergence
of $\tau$-GDA to DSE compared
to $\tau$-SGA.
However, in practice, it is numerically hard
to validate this theory in GAN.
In this section, we study
the training of orthogonal Wasserstein GAN
using stochastic $\tau$-GDA
and $\tau$-SGA.
We focus on stochastic gradients because
the number of training samples from a GAN generator
can be infinite.
We construct two non-standard examples with the following goals: (1).
illustrate a similar local convergence dilemma near DSE faced by $\tau$-GDA as $\tau$ varies
on both synthetic and real dataset.
(2). show how $\tau$-SGA may converge faster due to a smaller $\tau$.
(3). propose simple and clean benchmarks such that
theoretical study of Riemannian non-convex non-concave games
from GAN could be further developed in future works.

### 4.1 Setup of orthogonal Wasserstein GAN

In Wasserstein GAN,
we are interested in the following min-max problem

$$ $f(x,y)=\mathbb{E}(D_{y}(\phi_{data}))-\mathbb{E}(D_{y}(\phi_{x})).$ $$

The idea of GAN is to approximate the probability distribution
of $\phi_{data}\in\mathbb{R}^{d}$ using a generator $G_{x}$ parameterized by $x\in\mathcal{M}_{1}$.
The parameter $x$ is optimized so that
the distribution of $\phi_{x}$ approximates that of $\phi_{data}$
through the feedback of a discriminator $D_{y}:\mathbb{R}^{d}\to\mathbb{R}$,
parameterized by $y\in\mathcal{M}_{2}$.
We assume that $\mathcal{M}_{1}$ is Euclidean and $\mathcal{M}_{2}$ is Riemannian.
The manifold constraint on $y$ can ensure that the discriminator family is restricted to a sub-class of 1-Lipschitz continuous functions (up to a constant scaling).
The expectations in $f$ are taken with respect to
the random variable $\phi_{data}$ and $\phi_{x}$.(^2^22In the stochastic-gradient setting, the expectation of $D_{y}(\phi_{data})$ (resp. $D_{y}(\phi_{x})$) is estimated at each iteration of an algorithm, using a batch of samples of data in the training set (resp. a batch of samples of $Z$). The number of samples in a batch is called batch size.)

We consider unconditional GANs where a generative model $\phi_{x}=G_{x}(Z)$
transforms a random noise $Z\in\mathbb{R}^{p}$ to $\phi_{x}\in\mathbb{R}^{d}$.
We study two examples of $\phi_{data}$
using low-dimensional generators (dimension $p<d$).
In the second example, the discriminators is build upon the wavelet scattering transform
to capture discriminative information in natural images (Bruna & Mallat, ).
More details are given in Appendix L.

To study the local convergence of stochastic $\tau$-GDA and $\tau$-GDA,
we first obtain a reasonable solution of $(x,y)$ in terms of model quality using existing methods
(see Appendix L.5). We then evaluate these methods initialized from this solution.

##### Gaussian distribution by linear generator

We want to model $\phi_{data}\sim\mathcal{N}(0,\Sigma)$ on dimension $d=5$ using a PCA-like model with dimension $p=4$. The covariance matrix $\Sigma$ is diagonal with a small eigenvalue, i.e. $(1,2^{2},3^{2},4^{2},0.01)$. We choose
the generator $G_{x}(Z)=A_{x}Z$, where $A_{x}\in\mathbb{R}^{d\times p}$ and $Z$ is isotropic normal.
For the discriminator,
we use the Stiefel manifold $St(k,d)$ of $k$ orthogonal vectors on $\mathbb{R}^{d}$ to construct

$$ $D_{y}(\phi)=\langle v_{y},\sigma(W_{y}\phi)\rangle,$ $$

where $W_{y}\in St(k,d)$, $v_{y}\in St(1,k)$ and
$\sigma$ is a smooth non-linearity. We set $k=d=5$.
As $A_{x}$ is the only generator parameter, $x:=(A_{x})\in\mathcal{M}_{1}=\mathbb{R}^{d\times p}$. The discriminator parameter $y:=(W_{y},v_{y})\in\mathcal{M}_{2}=St(k,d)\times St(1,k)$.

##### MNIST (digit 0) by DCGAN generator

We aim to model the digit 0 from the MNIST
dataset ($d=28\times 28$)
using a DCGAN generator ($p=128$) in WGAN-GP (Gulrajani et al., ).
The space $\mathcal{M}_{1}$ is Euclidean with dimension $d_{1}=1556673$.
Unlike the CNN discriminator in WGAN-GP,
we consider a hybrid Scattering CNN discriminator (Oyallon et al., )
which has only one trainable layer,

$$ $D_{y}(\phi)=\langle v_{y},\sigma(w_{y}\star P(\phi)+b_{y})\rangle.$ $$

The scattering transform $P:\mathbb{R}^{d}\to\mathbb{R}^{I\times n\times n}$ is 1-Lipschitz-continuous (Bruna & Mallat, ) and
it has no trainable parameter. The scattering feature $P(\phi)\in\mathbb{R}^{I\times n\times n}$ is computed from an image $\phi$ of size $\sqrt{d}\times\sqrt{d}$.
We then apply an orthogonal convolutional layer to $P(\phi)$ (Cisse et al., ).
It has the kernel orthogonality (Achour et al., ), by reshaping
$w_{y}$ (size $J\times I\times k\times k$) into a matrix $W_{y}\in\mathbb{R}^{J\times(Ik^{2})}$ such that $W_{y}\in St(J,Ik^{2})$. The bias parameter $b_{y}\in\mathbb{R}^{J}$ and the output of the convolutional layer
is reshaped into a vector in $\mathbb{R}^{JN^{2}}$.
As in the Gaussian case, we further assume $v_{y}\in St(1,JN^{2})$. In summary, $y:=(W_{y},b_{y},v_{y})\in\mathcal{M}_{2}=St(J,Ik^{2})\times\mathbb{R}^{J}
\times St(1,JN^{2})$.

### 4.2 Results on Gaussian distribution by linear generator

In Figure 2(a),
we show how the function value $f(x(t),y(t))$
changes over iteration.
It is estimated using a batch of $1000$ training samples.
To use a fixed learning rate for the discriminator,
we set $\gamma=0.02/\tau$ for each choice of $\tau\in\{1,10,100\}$.
We observe that when $\tau=1$,
$f(x(t),y(t))$ has a huge oscillation,
but when $\tau=10$ and $100$, it has a much smaller oscillation amplitude.
To investigate the underlying reason,
we compute in Figure 2(b),
the evolution of an “angle” quantity every one thousand iterations.
The angle at $(x,y)$ is defined as
$\frac{\langle v_{y},\delta(x,y)\rangle}{\|\delta(x,y)\|}$,
where $\delta(x,y)=\mathbb{E}(\sigma(W_{y}\phi_{data}))-\mathbb{E}(\sigma(W_{y}\phi_{
x}))$ is estimated from the same batch of samples during the training.
When the angle is close to $1$, it indicates that $v_{y}$ is solved to be optimal.
We observe that when $\tau=1$, the angle oscillates between positive
and negative values, suggesting that $v_{y(t)}$ is detached from
$\delta(x(t),y(t))$. Therefore the minimization of $f$ over $x$
does not get a good feedback to improve the model. This is
confirmed in Figure 2(c),
where we compute the EMD distance (Rubner et al., ) every one thousand iterations,
between the empirical measure of $\phi_{data}$ and $\phi_{x(t)}$
(using $2000$ validation samples).
We observe that the EMD distance increases over $t$ when $\tau=1$.
On the contrary, it stays close to one
when $\tau=10$ and $100$. When $\tau=100$,
we find that the covariance error
$\|\Sigma-A_{x(t)}A_{x(t)}^{\intercal}\|$
stays around $0.3$ over all iterations.
This error is between
the smallest eigenvalue of $\Sigma$
(which is $0.01$, the PCA error)
and the second smallest eigenvalue of $\Sigma$.
Therefore to use a large enough $\tau$ is crucial
to reduce a large model error,
as it can ensure a positive angle to measure a meaningful distance between $\phi_{x}$ and $\phi_{data}$.
This phenomenon is consistent with the local convergence
of $\tau$-GDA to DSE in Figure 1(a).

Figure: (a) $f(x(t),y(t))$
Refer to caption: extracted/5893895/plots/plot_results_tau_effect_gaussian5bis_f.png

**Table 2: Last iteration measures of stochastic $\tau$-GDA and $\tau$-SGA on the Wasserstein GAN of MNIST (digit 0). We report the $f$, angle, and FID scores computed at $(x(T),y(T))$. $\tau$-GDA is trained for $T=2\times 10^{5}$ iterations with $\gamma=0.05/\tau$. $\tau$-SGA is trained longer with $\gamma=0.01$, $\tau=5$, $\theta=0.075$. Only significant digits are reported with respect to the standard deviation (shown in $(\pm)$) of the statistical estimation of these quantities (see footnote 3 and 4).**
|  | $\tau$-GDA |  |  |  |
| --- | --- | --- | --- | --- |
| $\tau$ | $f$ | angle | FID (train) | FID (val) |
| $5$ | -0.09 | -0.26 | 21 | 29 ($\pm 0.9$) |
| $10$ | 0.09 | 0.11 | 34 | 36 ($\pm 0.6$) |
| $30$ | 0.036 | 0.4 | 13 | 20 ($\pm 0.95$) |

### 4.3 Results on MNIST (digit 0) by DCGAN generator

We apply the insights from the theoretical analysis
to improve the convergence of both $\tau$-GDA
and $\tau$-SGA, in order to obtain a good generative model $\phi_{x}$.
To measure its quality,
we compute the FID scores (Seitzer, ).(^3^33 It is compared with the training, validation and test samples of a dataset. To compute FID (train) (resp. FID (val)), we generate the same number of fake samples as the training samples (resp. as the validation samples), and report an average value by re-sampling the fake samples $10$ times.)
The results of $\tau$-GDA
and $\tau$-SGA using batch size $128$
are reported in Table 2.
For $\tau$-GDA, we vary $\tau\in\{5,10,30\}$.
We observe that when $\tau=30$,
$\tau$-GDA has more stable
dynamics than $\tau=5$ and $\tau=10$
because the angle
stays around a positive constant.(^4^44To compute the angle, we compute $\delta(x,y)=\mathbb{E}(\sigma(w_{y}\star P(\phi_{data})+b_{y}))-\mathbb{E}( \sigma(w_{y}\star P(\phi_{x})+b_{y}))$. We estimate the $f$ and angle values by re-sampling $\phi_{x}$ ten times (similar to the estimation of FID). When $\tau=5$ and $\tau=10$, sometimes we observe a positive angle at $t=T$ in $\tau$-GDA because it oscillates around zero over $t$.)
On the other hand, a larger $\tau$ tends to slowdown the reduction of $f$
as in Figure 1(a).
Furthermore, the FID scores at $T=3\times 10^{6}$
are only slightly improved compared to those at $T=2\times 10^{5}$. Facing such a dilemma,
we evaluate the performance of $\tau$-SGA using
$\tau=5$ such that $\tau\gamma$ remains the same.
This is the case where $\tau$-GDA is not convergent due to
oscillating angles.
We find that both the angle and the FID scores are significantly improved
with a suitable choice of $\theta$ and $T$ in $\tau$-SGA.
We also perform a similar study on a subset of the Fashion-MNIST dataset,
which shows a consistent result (see Appendix L.6).
In terms of the computational time,
we find that $\tau$-GDA and $\tau$-SGA can
reach similar FID (test) scores after a same amount of training time (see Appendix K.2).

## 5 Conclusion

In this article, we provide sufficient conditions
for the local convergence
of $\tau$-GDA and $\tau$-SGA to differential equilibrium
on Riemannian manifold, using a classical fixed point theorem
based on spectral analysis.
This method allows one to reduce the problem into Euclidean space
using a local coordinate chart.
We obtain a linear local convergence of $\tau$-GDA to DSE and DNE
whose rate depends on the spectral radius of the Jacobian matrix of its update rule computed at the equilibrium.
To improve this convergence rate, $\tau$-SGA is developed and
extended to Riemannian manifold for the first time.
Based on the asymptotic analysis when $\tau$ is large,
we find that sometimes an improved rate of $\tau$-SGA to DSE can indeed be obtained.
The methodology could be served as a basis to analyze
other simultaneous algorithms such as the Riemannian
Hamiltonian method (Han et al., ).

To show the relevance of our results to GAN,
we study the behavior of stochastic $\tau$-GDA and $\tau$-SGA
on simple benchmarks of orthogonal Wasserstein GANs.
Even though our current theory does not apply to analyze
these stochastic algorithms, we observe a consistent behavior similar to
their deterministic convergence to DSE.
This suggests that DSE seems to be
a suitable solution set to which
existing min-max algorithms of GAN converge,
as observed in other games (non-zero sum) arising from
NS-GAN and WGAN-GP (Berard et al., ).

## Appendix A Definition of DSE

Let
$\mathcal{X}(\mathcal{M}_{1})$ (resp. $\mathcal{X}(\mathcal{M}_{2})$)
be the set of continuously differentiable vector fields on $\mathcal{M}_{1}$ (resp. $\mathcal{M}_{2}$).
Since $f$ is twice continuously differentiable,
we have for any $y\in\mathcal{M}_{2}$ (resp. $x\in\mathcal{M}_{1}$),
$\mbox{grad}_{x}f(\cdot,y)\in\mathcal{X}(\mathcal{M}_{1})$
(resp. $\mbox{grad}_{y}f(x,\cdot)\in\mathcal{X}(\mathcal{M}_{2})$).
We write $g_{1}$ and $g_{2}$ to denote
the Riemannian metric on $\mathcal{M}_{1}$ and $\mathcal{M}_{2}$.

The main idea of the proof is to use the local coordinate chart
$(O_{1}\times O_{2},\varphi_{1}\times\varphi_{2})$ to represent the smooth vector fields
$(x,y)\mapsto\mbox{grad}_{x}(x,y)$ and $(x,y)\mapsto\mbox{grad}_{y}(x,y)$ around the point $(x^{\ast},y^{\ast})$, so as to show the equivalence
between (4)-(6)
and (1)-(3).

For each $x\in O_{1}$ (resp. $y\in O_{2}$), let
$\{E_{1,i}(x)\}_{i\leq d_{1}}$ (resp. $\{E_{2,j}(y)\}_{j\leq d_{2}}$)
be the canonical basis of the tangent space $T_{x}\mathcal{M}_{1}$ (resp. $T_{y}\mathcal{M}_{2}$),
defined by the tangent map
$D\varphi^{-1}_{1}(\varphi_{1}(x))[e_{1,i}]$ of the canonical basis $\{e_{1,i}\}_{i\leq d_{1}}$ on $\mathbb{R}^{d_{1}}$
(resp.
$D\varphi^{-1}_{2}(\varphi_{2}(y))[e_{2,j}]$ of $\{e_{2,j}\}_{j\leq d_{2}}$ on $\mathbb{R}^{d_{2}}$).

Let the local coordinate of $(x,y)\in O_{1}\times O_{2}$
be $(u_{1},u_{2})=(\varphi_{1}(x),\varphi_{2}(y))\in\bar{O}_{1}\times\bar{O}_{2}=
\varphi_{1}(O_{1})\times\varphi_{2}(O_{2})\subset\mathbb{R}^{d_{1}}\times
\mathbb{R}^{d_{2}}$.
We can represent

$$ $\displaystyle\xi_{1}(x,y)=\mbox{grad}_{x}f(x,y)$ $\displaystyle=\sum_{i\leq d_{1}}\xi_{1,i}(x,y)E_{1,i}(x),$ (15) $\displaystyle\xi_{2}(x,y)=\mbox{grad}_{y}f(x,y)$ $\displaystyle=\sum_{j\leq d_{2}}\xi_{2,j}(x,y)E_{2,j}(y).$ (16) $$

From (15) and (16),
the local coordinates of the Riemannian gradients on $\bar{O}_{1}\times\bar{O}_{2}$ are

$$ $\displaystyle\bar{\xi}_{1}(u_{1},u_{2})$ $\displaystyle=\sum_{i\leq d_{1}}\xi_{1,i}(\varphi^{-1}_{1}(u_{1}),\varphi^{-1} _{2}(u_{2}))e_{1,i},$ (17) $\displaystyle\bar{\xi}_{2}(u_{1},u_{2})$ $\displaystyle=\sum_{j\leq d_{2}}\xi_{2,j}(\varphi^{-1}_{1}(u_{1}),\varphi^{-1} _{2}(u_{2}))e_{2,j}.$ (18) $$

##### Equivalence between ( 4 ) and ( 1 )

Let the local Riemannian metric matrix at $x=\varphi^{-1}_{1}(u_{1})\in O_{1}$
and $y=\varphi^{-1}_{2}(u_{2})\in O_{2}$ be

$$ $\bar{g}_{1}(u_{1})=(g_{1}(E_{1,i}(x),E_{1,i^{\prime}}(x)))_{i,i^{\prime}\leq d _{1}},\quad\bar{g}_{2}(u_{2})=(g_{2}(E_{2,j}(y),E_{2,j^{\prime}}(y)))_{j,j^{ \prime}\leq d_{2}}.$ (19) $$

From (Absil et al., , chap. 3.6), we have for $(u_{1},u_{2})\in\bar{O}_{1}\times\bar{O}_{2}$,

$$ $\displaystyle\bar{\xi}_{1}(u_{1},u_{2})=\bar{g}_{1}(u_{1})^{-1}\cdot\partial_{ u_{1}}\bar{f}(u_{1},u_{2}),$ (20) $\displaystyle\bar{\xi}_{2}(u_{1},u_{2})=\bar{g}_{2}(u_{2})^{-1}\cdot\partial_{ u_{2}}\bar{f}(u_{1},u_{2}).$ (21) $$

The equivalence between the DSE condition (4) and (1) follows from (20) and (21), as $\bar{\xi}_{1}(u_{1}^{\ast},u_{2}^{\ast})=0$ i.f.f. (if and only if) $\mbox{grad}_{x}f(x^{\ast},y^{\ast})=0$ (similarly for the relationship between $\bar{\xi}_{2}$ and $\mbox{grad}_{y}f$).

##### Equivalence between ( 5 ) and ( 2 )

Let $\eta^{\ast}\in T_{y^{\ast}}\mathcal{M}_{2}$,
with its local coordinate $\bar{\eta}^{\ast}=D\varphi_{2}(y^{\ast})[\eta^{\ast}]$.
To show that the positive definiteness of
$-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$
is equivalent to
the positive definiteness
of $-\partial^{2}_{u_{2}u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})$,
it is sufficient to verify that

$$ $\langle\bar{\eta}^{\ast},\partial^{2}_{u_{2}u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{ \ast})\bar{\eta}^{\ast}\rangle=\langle\eta^{\ast},\mbox{Hess}_{y}f(x^{\ast},y^ {\ast})[\eta^{\ast}]\rangle_{y^{\ast}},$ (22) $$

where $\langle\cdot,\cdot\rangle$ denotes the inner product on Euclidean space,
and $\langle\cdot,\cdot\rangle_{y^{\ast}}$ denotes the Riemannian inner product on $T_{y^{\ast}}M$.

Recall that the Riemannian Hessian on $\mathcal{M}_{2}$ is defined
by the Riemannian connection $\nabla_{y}:\mathcal{X}(\mathcal{M}_{2})\times\mathcal{X}(\mathcal{M}_{2})\to
\mathcal{X}(\mathcal{M}_{2})$
with $\mbox{Hess}_{y}f(x,y)=\nabla_{y}\xi_{2}(x,y)$. By the R-$linearity$ and Leibniz law of the connection, we have

$$ $\displaystyle\nabla_{y}\xi_{2}(x^{\ast},y^{\ast})[\eta^{\ast}]$ $\displaystyle=\nabla_{y}(\sum_{j}\xi_{2,j}(x^{\ast},y^{\ast})E_{2,j}(y^{\ast}) )[\eta^{\ast}]$ $\displaystyle=\sum_{j}\nabla_{y}(\xi_{2,j}(x^{\ast},y^{\ast})E_{2,j}(y^{\ast}) )[\eta^{\ast}]$ $\displaystyle=\sum_{j}(D_{y}\xi_{2,j}(x^{\ast},y^{\ast})[\eta^{\ast}])E_{2,j}( y^{\ast})+\xi_{2,j}(x^{\ast},y^{\ast})\nabla_{y}(E_{2,j}(y^{\ast}))[\eta^{\ast}]$ (23) $\displaystyle=\sum_{j}\left(\frac{\partial{\bar{\xi}_{2,j}}}{\partial u_{2}}(u _{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}\right)E_{2,j}(y^{\ast}).$ (24) $$

We have obtained (24) from (23) because any DSE $(x^{\ast},y^{\ast})$ is a critical point of $f$, therefore $\xi_{2,j}(x^{\ast},y^{\ast})=0$ for each $j\leq d_{2}$. Moreover, the tangent map
$D_{y}\xi_{2,j}(x^{\ast},y^{\ast})[\eta^{\ast}]=\frac{\partial{\bar{\xi}_{2,j}}
}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}$ by definition (Absil et al., , chap. 3.5.4).

Therefore the local coordinate of the tangent vector $\nabla_{y}\xi_{2}(x^{\ast},y^{\ast})[\eta^{\ast}]\in T_{y^{\ast}}\mathcal{M}_{2}$ is $\frac{\partial{\bar{\xi}_{2}}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{
\eta}^{\ast}\in\mathbb{R}^{d_{2}}$.
It follows from (21),(24),
and the Riemannian inner product (Absil et al., , chap. 3.6)
that

$$ $\displaystyle\langle\eta^{\ast},\mbox{Hess}_{y}f(x^{\ast},y^{\ast})[\eta^{\ast }]\rangle_{y^{\ast}}$ $\displaystyle=\left\langle\bar{\eta}^{\ast},\bar{g}_{2}(u_{2}^{\ast})\cdot \frac{\partial{\bar{\xi}_{2}}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{ \eta}^{\ast}\right\rangle$ (25) $\displaystyle=\left\langle\bar{\eta}^{\ast},\partial^{2}_{u_{2}u_{2}}\bar{f}(u _{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}\right\rangle.$ $$

Therefore (22) holds.
This verifies the equivalence between (5) and (2).

##### Equivalence between ( 6 ) and ( 3 )

Assume $(x^{\ast},y^{\ast})$ is a DSE of $f$,
from the equivalence between
(4),(5) and (1),(2),
we deduce that $(x^{\ast},y^{\ast})$ (resp. $(u_{1}^{\ast},u_{2}^{\ast})$) is a critical
point of $f$ (resp. $\bar{f}$).
Moreover, $\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ and $\partial^{2}_{u_{2}u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})$
are negative definite.

To simplify the notation, let

$$ $\displaystyle A=-\partial^{2}_{u_{2}u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast}), \quad B=\partial^{2}_{u_{2}u_{1}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast}),\quad C= \partial^{2}_{u_{1}u_{1}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast}).$ (26) $$

The condition (6) is equivalent to the positive definiteness
of the matrix $C+BA^{-1}B^{\intercal}$.

Let $\delta^{\ast}\in T_{x^{\ast}}\mathcal{M}_{1}$,
with its local coordinate $\bar{\delta}^{\ast}=D\varphi_{1}(x^{\ast})[\delta^{\ast}]$.
As in (22), it suffices to verify that

$$ $\langle\bar{\delta}^{\ast},(C+BA^{-1}B^{\intercal})\bar{\delta}^{\ast}\rangle= \langle\delta^{\ast},[\mbox{Hess}_{x}f-\mbox{grad}^{2}_{yx}f\cdot(\mbox{Hess}_ {y}f)^{-1}\cdot\mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast})[\delta^{\ast}]\rangle _{x^{\ast}}.$ $$

Based on the proof of (22),
we also have $\langle\bar{\delta}^{\ast},C\bar{\delta}^{\ast}\rangle=\langle\delta^{\ast},
\mbox{Hess}_{x}f(x^{\ast},y^{\ast})[\delta^{\ast}]\rangle_{x^{\ast}}$. It remains to verify that

$$ $\langle\bar{\delta}^{\ast},BA^{-1}B^{\intercal}\bar{\delta}^{\ast}\rangle= \langle\delta^{\ast},[-\mbox{grad}^{2}_{yx}f\cdot(\mbox{Hess}_{y}f)^{-1}\cdot \mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast})[\delta^{\ast}]\rangle_{x^{\ast}}.$ (27) $$

Let $\eta^{\ast}=[\mbox{Hess}_{y}f^{-1}\cdot\mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast
})[\delta^{\ast}]$, we compute the local coordinate $\bar{\eta}^{\ast}$ of $\eta^{\ast}$ by converting
the following equation

$$ $\mbox{Hess}_{y}f(x^{\ast},y^{\ast})[\eta^{\ast}]=\mbox{grad}^{2}_{xy}f(x^{\ast },y^{\ast})[\delta^{\ast}],$ (28) $$

into local coordinate. From (24), we have
$\mbox{Hess}_{y}f(x^{\ast},y^{\ast})[\eta^{\ast}]=\sum_{j}\left(\frac{\partial{
\bar{\xi}_{2,j}}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}
\right)E_{2,j}(y^{\ast})$.

By (16) and (18), we obtain

$$ $\displaystyle\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]$ $\displaystyle=D_{x}\mbox{grad}_{y}f(x^{\ast},y^{\ast})[\delta^{\ast}]=\sum_{j} D_{x}\xi_{2,j}(x^{\ast},y^{\ast})[\delta^{\ast}]E_{2,j}(y^{\ast})$ $\displaystyle=\sum_{j}\left(\frac{\partial\bar{\xi}_{2,j}}{\partial u_{1}}(u_{ 1}^{\ast},u_{2}^{\ast})\bar{\delta}^{\ast}\right)E_{2,j}(y^{\ast}).$ (29) $$

It follows that (28) is equivalent to

$$ $\displaystyle\frac{\partial{\bar{\xi}_{2}}}{\partial u_{2}}(u_{1}^{\ast},u_{2} ^{\ast})\bar{\eta}^{\ast}=\frac{\partial\bar{\xi}_{2}}{\partial u_{1}}(u_{1}^{ \ast},u_{2}^{\ast})\bar{\delta}^{\ast}.$ (30) $$

We can compute the right hand side of (27) in the local coordinate as in (29) and (25),

$$ $\displaystyle\langle\delta^{\ast},-\mbox{grad}^{2}_{yx}f(x^{\ast},y^{\ast})[ \eta^{\ast}]\rangle_{x^{\ast}}$ $\displaystyle=-\langle\delta^{\ast},D_{y}\mbox{grad}_{x}f(x^{\ast},y^{\ast})[ \eta^{\ast}]\rangle_{x^{\ast}}$ $\displaystyle=-\left\langle\delta^{\ast},\sum_{i}\left(\frac{\partial\bar{\xi} _{1,i}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}\right)E_{1 ,i}(x^{\ast})\right\rangle_{x^{\ast}}$ $\displaystyle=-\left\langle\bar{\delta}^{\ast},\bar{g}_{1}(u_{1}^{\ast})\cdot \frac{\partial\bar{\xi}_{1}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{ \eta}^{\ast}\right\rangle.$ $$

The left hand side of (27)
can computed based on (20) and (21),
which results in

$$ $\displaystyle A$ $\displaystyle=-\bar{g}_{2}(u_{2}^{\ast})\cdot\frac{\partial{\bar{\xi}_{2}}}{ \partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast}),$ (31) $\displaystyle B$ $\displaystyle=\bar{g}_{1}(u_{1}^{\ast})\cdot\frac{\partial{\bar{\xi}_{1}}}{ \partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast}),$ (32) $\displaystyle C$ $\displaystyle=\bar{g}_{1}(u_{1}^{\ast})\cdot\frac{\partial{\bar{\xi}_{1}}}{ \partial u_{1}}(u_{1}^{\ast},u_{1}^{\ast}).$ $$

Furthermore, the symmetry of the Hessian matrix of $\bar{f}$ implies that

$$ $B^{\intercal}=\partial^{2}_{u_{1}u_{2}}\bar{f}(u_{1}^{\ast},u_{2}^{\ast})=\bar {g}_{2}(u_{2}^{\ast})\cdot\frac{\partial\bar{\xi}_{2}}{\partial u_{1}}(u_{1}^{ \ast},u_{2}^{\ast}).$ (33) $$

As a consequence, (30) implies that $\bar{\eta}^{\ast}=-A^{-1}B^{\intercal}\bar{\delta}^{\ast}$, and

$$ $\langle\bar{\delta}^{\ast},BA^{-1}B^{\intercal}\bar{\delta}^{\ast}\rangle=- \left\langle\bar{\delta}^{\ast},\bar{g}_{1}(u_{1}^{\ast})\cdot\frac{\partial{ \bar{\xi}_{1}}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast} \right\rangle.$ $$

Therefore (27) holds.

## Appendix B Implicit function theorem on Riemannian manifold

We use an implicit function theorem
which allows one to understand why $(x^{\ast},y^{\ast})$ is a local minimax point
in the manifold case.
This theorem implies the existence of a solution $(x,h(x))$ sufficiently
close to $(x^{\ast},y^{\ast})$ s.t. $\mbox{grad}_{y}f(x,h(x))=0$.
Moreover, due to the continuity of $\mbox{Hess}_{y}f$ on $\mathcal{M}_{1}\times\mathcal{M}_{2}$ and the continuity of $h$ near $x^{\ast}$,
$-\mbox{Hess}_{y}f(x,h(x))$ is positive definite,
provided that $x$ is close enough to $x^{\ast}$.
As a consequence, $h(x)$ is a unique strict local maximum of
$y\mapsto f(x,y)$ in this neighbor of DSE.

Denote the set of continuously differentiable vector fields on $\mathcal{M}_{2}$ by $\mathcal{X}(\mathcal{M}_{2})$.
Consider a parameterized vector field $\xi$ defined at each $x\in\mathcal{M}_{1}$ with
$\xi(x,\cdot)\in\mathcal{X}(\mathcal{M}_{2})$. Recall that the Riemannian connection on $\mathcal{M}_{2}$ is denoted by
$\nabla_{y}$.
By definition, for each $(x,y)\in\mathcal{M}_{1}\times\mathcal{M}_{2}$, $\xi(x,y)\in T_{y}\mathcal{M}_{2}$,
and $\nabla_{y}\xi(x,y)$ is a linear map from $T_{y}\mathcal{M}_{2}$ to $T_{y}\mathcal{M}_{2}$.

###### Theorem B.1 .

Assume $\xi$ is continuously differentiable on an open set $E\subset\mathcal{M}_{1}\times\mathcal{M}_{2}$.
Let $(x^{\ast},y^{\ast})\in E$ be a solution of

$$ $\xi(x,y)=0,\quad(x,y)\in E.$ $$

If $\nabla_{y}\xi(x^{\ast},y^{\ast})$ is invertible on $T_{y^{\ast}}\mathcal{M}_{2}$,
there exists an open set $U_{1}\times U_{2}\subset E$ and
an open set $W_{1}\subset U_{1}$ such that

$$ $\forall x\in W_{1},\exists!\,y\in U_{2}\quad\mbox{s.t.}\quad\xi(x,y)=0.$ $$

Let the unique $y=h(x)$, i.e. $h:W_{1}\to U_{2}$,
such that

$$ $y^{\ast}=h(x^{\ast}),\quad\xi(x,h(x))=0,\quad\forall x\in W_{1},$ $$

then $h$ is continuously differentiable on $W_{1}$.
Let $\delta\in T_{x^{\ast}}\mathcal{M}_{1}$,
and denote the tangent map of $h$ at $x^{\ast}$ by $D_{x}h(x^{\ast})$,
we have

$$ $D_{x}h(x^{\ast})[\delta]=-\nabla_{y}\xi(x^{\ast},y^{\ast})^{-1}\cdot D_{x}\xi( x^{\ast},y^{\ast})[\delta].$ $$

We may prove this result by using a local coordinate chart
around the point $(x^{\ast},y^{\ast})$ to represent
a smooth vector field, and by the proof technique of the implicit function theorem in Euclidean space (Rudin, , Theorem 9.28).

###### Theorem B.1 .

## Appendix C Proof of Proposition 2.1

We verify the intrinsic definition of DSE for $(x^{\ast},y^{\ast})$.
For the first-order condition in (1), we compute based on
Boumal (, Proposition 3.61),

$$ $\partial_{x}f(x,y)=A^{\intercal}y,\quad\mbox{grad}_{y}f(x,y)=(I-yy^{\intercal} )(Ax-b).$ (34) $$

From the identity involving the pseudo-inverse
$A^{\intercal}AA^{+}=A^{\intercal}$,
we deduce that $A^{\intercal}y^{\ast}=0$.
As $b\not\in\mbox{Range}(A)$, we have
$Ax^{\ast}-b\neq 0$ and therefore
$y^{\ast}$ is parallel to the (non-zero) vector $Ax^{\ast}-b$.
From above, the first-order condition holds, i.e.
$\partial_{x}f(x^{\ast},y^{\ast})=0,\mbox{grad}_{y}f(x^{\ast},y^{\ast})=0$.

From (34), we first verify the second-order condition (2).
From Boumal (, Corollary 5.16),
we compute for $\eta\in T_{y}\mathcal{M}_{2}$,

$$ $\displaystyle\mbox{Hess}_{y}f(x,y)[\eta]$ $\displaystyle=(I-yy^{\intercal})(\langle\partial_{y}\mbox{grad}_{y}f(x,y),\eta\rangle)$ $\displaystyle=(I-yy^{\intercal})(-(y\eta^{\intercal}+\eta y^{\intercal})(Ax-b))$ $\displaystyle=\langle y,Ax-b\rangle(-\eta).$ $$

We verify that $-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ is d.p, because for non-zero $\eta^{\ast}\in T_{y^{\ast}}\mathcal{M}_{2}$, we have
$\|\eta^{\ast}\|^{2}>0$. Moreover, $Ax^{\ast}-b\neq 0$, therefore

$$ $\langle-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})[\eta^{\ast}],\eta^{\ast}\rangle= \langle y^{\ast},Ax^{\ast}-b\rangle\|\eta^{\ast}\|^{2}=\|Ax^{\ast}-b\|\|\eta^{ \ast}\|^{2}>0.$ $$

We now check the second-order condition (3).
It is clear that $\partial^{2}_{xx}f(x,y)=0$.
We next show $\delta^{\ast}\mapsto\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]$ is an injection from $T_{x^{\ast}}\mathcal{M}_{1}$ to $T_{y^{\ast}}\mathcal{M}_{2}$. We compute from (34),

$$ $\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]=\sum_{i\leq d_{1}} \partial_{x_{i}}\mbox{grad}_{y}f(x^{\ast},y^{\ast})\delta^{\ast}_{i}=(I-y^{ \ast}y^{\ast\intercal})A\delta^{\ast}.$ $$

If $\delta^{\ast}\neq 0$, then $A\delta^{\ast}\neq 0$ because $\mbox{Ker}(A)=\{0\}$. But $A\delta^{\ast}$ is not parallel
to $y^{\ast}$, since $y^{\ast}$ is along the direction
$Ax^{\ast}-b$ which is not in the range of $A$. This proves that $\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]\neq 0$ if $\delta^{\ast}\neq 0$.

The above injection property implies that for $\delta^{\ast}\neq 0$, $\eta^{\ast}=\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]\neq 0$. As $-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ is d.p, we use the symmetry of the Riemannian cross-gradients (see Proposition 3.1) to obtain: If $\delta^{\ast}\neq 0$, then

$$ $\langle\delta^{\ast},[\mbox{Hess}_{x}f-\mbox{grad}^{2}_{yx}f\cdot(\mbox{Hess}_ {y}f)^{-1}\cdot\mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast})[\delta^{\ast}]\rangle =\langle\eta^{\ast},-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})^{-1}[\eta^{\ast}] \rangle>0.$ $$

Therefore (3) holds.

## Appendix D Proof of Proposition 2.2

As in the proof of Proposition 2.1 in Appendix C, we have

$$ $\partial_{x}f(x,y)=A^{\intercal}y-\kappa A^{\intercal}Ax,\quad\mbox{grad}_{y}f (x,y)=(I-yy^{\intercal})(Ax-b).$ (35) $$

To show the existence of DSE, we shall construct a solution
of the form $(x,y^{\ast}_{x})$, where
$\partial_{x}f(x,y^{\ast}_{x})=0$.
We consider $y^{\ast}_{x}=(Ax-b)/\|Ax-b\|$
so that $\mbox{grad}_{y}f(x,y^{\ast}_{x})=0$.
We verify that $\mbox{Ker}(A)=\{0\}$ implies that $A^{\intercal}A$ is p.d.
It follows that $A^{+}=(A^{\intercal}A)^{-1}A^{\intercal}$,
and that the condition $\partial_{x}f(x,y^{\ast}_{x})=0$ is equivalent to

$$ $\frac{A^{\intercal}(Ax-b)}{\|Ax-b\|}=\kappa A^{\intercal}Ax\quad \Leftrightarrow\quad(1-\kappa\|Ax-b\|)x=A^{+}b.$ $$

For $(x^{\ast},y^{\ast})=(x^{\ast},y^{\ast}_{x^{\ast}})$ to be a DSE,
we assume $x^{\ast}=cA^{+}b$
and we want to find a $c\in\mathbb{R}$ such that

$$ $F(\kappa,c)=c(1-\kappa\|cAA^{+}b-b\|)-1=0.$ $$

From Proposition 2.1, we have $F(0,1)=0$.
We next apply the implicit function theorem (c.f. Theorem B.1)
to show the existence of $c$ close to $1$ if $\kappa$ is close to 0 0.
It suffices to verify that $\frac{\partial F}{\partial c}(0,1)=1$ based on

$$ $\frac{\partial F}{\partial c}(\kappa,c)=\left(1-\kappa\|cAA^{+}b-b\|\right)+c \left(-\kappa\frac{\langle AA^{+}b,cAA^{+}b-b\rangle}{\|cAA^{+}b-b\|}\right).$ $$

This suggests that there is
an implicit function $h$ defined on an open neighbor $W_{1}$ of $\kappa=0$ ($0\in W_{1}$)
such that $F(\kappa,h(\kappa))=0$. Moreover, $h(0)=1$.
We next check that if $c=h(\kappa)$ for some range $0<\kappa<\kappa_{0}$,
then $(x^{\ast},y^{\ast})$ defined in the statement is a DSE of $f$.
Since (35) hold at $(x^{\ast},y^{\ast})$,
it suffices to verify the second-order condition of DSE.

We verify the second-order condition (2)
by using the same proof of Proposition 2.1.
To check the second-order condition (3),
we first compute from (35)

$$ $\eta^{\ast}=\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]=\sum_{i \leq d_{1}}\partial_{x_{i}}\mbox{grad}_{y}f(x^{\ast},y^{\ast})\delta^{\ast}_{i }=(I-y^{\ast}y^{\ast\intercal})A\delta^{\ast}.$ $$

By following the proof of Proposition 2.1,
we then compute

$$ $\displaystyle S(\kappa)$ $\displaystyle=\langle\delta^{\ast},[\mbox{Hess}_{x}f-\mbox{grad}^{2}_{yx}f \cdot(\mbox{Hess}_{y}f)^{-1}\cdot\mbox{grad}^{2}_{xy}f](x^{\ast},y^{\ast})[ \delta^{\ast}]\rangle$ $\displaystyle=-\kappa\|A\delta^{\ast}\|^{2}+\langle\eta^{\ast},-\mbox{Hess}_{y }f(x^{\ast},y^{\ast})^{-1}[\eta^{\ast}]\rangle$ $\displaystyle=-\kappa\|A\delta^{\ast}\|^{2}+\frac{1}{\|Ax^{\ast}-b\|}\|\eta^{ \ast}\|^{2}.$ $$

It is clear that $S(\kappa)$ is a continuous function of $\kappa$
because $x^{\ast}=h(\kappa)A^{+}b$ and
$y^{\ast}=(Ax^{\ast}-b)/\|Ax^{\ast}-b\|$ are continuous with respect to $\kappa$.
Moreover, $S(0)>0$ from the proof of Proposition 2.1
(due to the injection property of $\delta^{\ast}\mapsto\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})[\delta^{\ast}]$).
Therefore there exists $\kappa_{0}>0$ so that $S(\kappa)>0$ for $0<\kappa<\kappa_{0}$ (i.e. (3) holds).

## Appendix E Proof of Proposition 2.3

We verify the intrinsic definition of DNE for $(x^{\ast},y^{\ast})$.
For the first-order condition (7), we compute

$$ $\partial_{x}f(x,y)=A^{\intercal}(Ax+y-b),\quad\mbox{grad}_{y}f(x,y)=(I-yy^{ \intercal})(Ax+y-b).$ (36) $$

From the identity involving the pseudo-inverse
$A^{\intercal}AA^{+}=A^{\intercal}$,
we deduce that $A^{\intercal}(Ax^{\ast}+y^{\ast}-b)=0$.
As $b\not\in\mbox{Range}(A)$, we have
$Ax^{\ast}-b\neq 0$ and therefore
$y^{\ast}$ is parallel to the vector $Ax^{\ast}-b$.
From above, the first-order condition in (1) holds.

From (36), we next verify the second-order condition (8).
It is clear that $\partial^{2}_{xx}f(x,y)=A^{\intercal}A$ is p.d. From Boumal (, Corollary 5.16),
we compute for $\eta\in T_{y}\mathcal{M}_{2}$,

$$ $\displaystyle\mbox{Hess}_{y}f(x,y)[\eta]$ $\displaystyle=(I-yy^{\intercal})(\langle\partial_{y}\mbox{grad}_{y}f(x,y),\eta\rangle)$ $\displaystyle=(I-yy^{\intercal})(\eta-(y\eta^{\intercal}+\eta y^{\intercal})( Ax+y-b))$ $\displaystyle=\eta(1-\langle y,Ax+y-b\rangle)$ $\displaystyle=\langle y,Ax-b\rangle(-\eta).$ $$

We verify that $-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$ is d.p, because for non-zero $\eta^{\ast}\in T_{y^{\ast}}\mathcal{M}_{2}$, we have
$\|\eta^{\ast}\|^{2}>0$. Moreover, $Ax^{\ast}-b\neq 0$, therefore

$$ $\langle-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})[\eta^{\ast}],\eta^{\ast}\rangle= \langle y^{\ast},Ax^{\ast}-b\rangle\|\eta^{\ast}\|^{2}=\|Ax^{\ast}-b\|\|\eta^{ \ast}\|^{2}>0.$ $$

Therefore (8) holds.

## Appendix F Proof of Theorem 3.1

The proof is illustrated in Figure 3,
which contains three main steps:
(1) identify a set $\bar{S}$ where
the local the local update rule $\bar{\mathbf{T}}$
is well-defined.
(2) adapt the proof technique
of the classical Ostrowski Theorem
to identity a local stable region $\bar{S}_{\delta}$
and pull it back to the manifold (i.e. $S_{\delta}$).
(3) relate the Euclidean distance in the local coordinate
to the Riemannian distance on $S_{\delta}$ to
establish the linear convergence rate.

Figure: (a)

##### Preliminary

Let $\varphi=\varphi_{1}\times\varphi_{2}$.
Let us consider the set $S=\mathbf{T}^{-1}(O_{1}\times O_{2})\cap(O_{1}\times O_{2})$
and its local coordinate domain $\bar{S}=(\varphi_{1}\times\varphi_{2})(S)\subset\bar{O}_{1}\times\bar{O}_{2}$.
By definition, $(x,y)\in S\subset O_{1}\times O_{2}$ and $\mathbf{T}(x,y)\in O_{1}\times O_{2}$.
Therefore the induced dynamics $\bar{\mathbf{T}}$ is well-defined on $\bar{S}$, i.e. $\forall(u_{1},u_{2})\in\bar{S}$,

$$ $\bar{\mathbf{T}}(u_{1},u_{2})=\varphi(\mathbf{T}(\varphi^{-1}_{1}(u_{1}), \varphi^{-1}_{2}(u_{2}))).$ $$

Note that $S$ and $\bar{S}$ are non-empty open sets
since $\mathbf{T}$ is continuous (by the continuity of
the vector fields $\xi_{1}$, $\xi_{2}$ and the retractions
$\mathcal{R}_{1}$ and $\mathcal{R}_{2}$) and $(x^{\ast},y^{\ast})$ is a fixed point of $\mathbf{T}$.

##### Local stable region

We aim to identify a local stable region $\bar{S}_{\delta}\subset\bar{S}$
around $(u_{1}^{\ast},u_{2}^{\ast})$ (the local coordinate of $(x^{\ast},y^{\ast})$),
so that if $(u_{1},u_{2})\in\bar{S}_{\delta}$, then
$\bar{\mathbf{T}}(u_{1},u_{2})\in\bar{S}_{\delta}\subset\bar{O}_{1}\times\bar{O
}_{2}$.
Let $S_{\delta}=\varphi^{-1}(\bar{S}_{\delta})$,
and assume $(x(0),y(0))\in S_{\delta}$,
then by recursion $\bar{\mathbf{T}}(\varphi_{1}(x(t)),\varphi_{2}(y(t)))$ is always well defined
since the sequence $(x(t),y(t))\in S_{\delta}\subset S,\forall t\geq 1$.

To construct such a region, the key is to verify that the spectral radius $\rho$
of the Jacobian matrix $\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$
is strictly smaller than one. In Appendix F.1,
we verify that the eigenvalues of $\mathbf{T}^{\prime}(x^{\ast},y^{\ast})$
are the same as $\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$.
Therefore according to our assumption,
$\rho:=\rho(\mathbf{T}^{\prime}(x^{\ast},y^{\ast})<1$.
From the proof of Ostrowski Theorem,
for an arbitrary $\epsilon>0$,
there exists a norm $\|\cdot\|_{\epsilon}$ on $\mathbb{R}^{d_{1}}\times\mathbb{R}^{d_{2}}$,
such that

$$ $\|\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})\|_{\epsilon}\leq\rho+ \epsilon/2.$ $$

Furthermore, due to the differentiability of $\bar{\mathbf{T}}$ at the fixed point,
there exists $\delta>0$, and
$\bar{S}_{\delta}^{\prime}=\{(u_{1},u_{2})|\|(u_{1},u_{2})-(u_{1}^{\ast},u_{2}^
{\ast})\|_{\epsilon}<\delta\}\subset\bar{S}$, such that

$$ $\displaystyle\|\bar{\mathbf{T}}(u_{1},u_{2})-(u_{1}^{\ast},u_{2}^{\ast})\|_{\epsilon}$ $\displaystyle\leq(\|\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})\|_{ \epsilon}+\epsilon/2)\|(u_{1},u_{2})-(u_{1}^{\ast},u_{2}^{\ast})\|_{\epsilon}$ $\displaystyle\leq(\rho+\epsilon)\|(u_{1},u_{2})-(u_{1}^{\ast},u_{2}^{\ast})\|_ {\epsilon},\quad\forall(u_{1},u_{2})\in\bar{S}_{\delta}^{\prime}.$ (37) $$

As the set $\bar{S}_{\delta}^{\prime}$ is open,
we can identify an open geodesically convex subset of
$\varphi^{-1}(\bar{S}_{\delta}^{\prime})$
as the local stable region $S_{\delta}$.
The set $S_{\delta}$ is also homeomorphic to a ball, i.e. without any hole.

##### Locally convergent with rate ρ 𝜌 \rho italic_ρ

Let’s choose $\epsilon$ so that $\rho+\epsilon<1$,
then from (F),
if $(u_{1}(0),u_{2}(0))\in\bar{S}_{\delta}$, the sequence $(u_{1}(t),u_{2}(t))$ stays in $\bar{S}_{\delta}$ and
converges to $(u_{1}^{\ast},u_{2}^{\ast})$ as $t\to\infty$.
As $\varphi=\varphi_{1}\times\varphi_{2}$ is a continuous bijection (homeomorphism)
from $O_{1}\times O_{2}$ to $\bar{O}_{1}\times\bar{O}_{2}$, we have equivalently that
if $(x(0),y(0))\in S_{\delta}=\varphi^{-1}(\bar{S}_{\delta})$,
the sequence $(x(t),y(t))$ will stay in $S_{\delta}$ and
converges to $(x^{\ast},y^{\ast})$.

The linear convergence rate in the local coordinate
system in (F) can be used
to control the Riemannian distance $d(t)$
between $\mathbf{p}(t)=(x(t),y(t))\in\mathcal{M}=\mathcal{M}_{1}\times\mathcal{M}_{2}$
and $\mathbf{p}^{\ast}=(x^{\ast},y^{\ast})\in\mathcal{M}$.
As $S_{\delta}$ is geodesically convex,
the Riemannian distance equals to the geodesic distance.
Let $\Gamma_{\delta}(t)$ be the set of piece-wise smooth curves
restricted on $S_{\delta}$ with initial point $\mathbf{p}(t)$ at time 0 0 and
last point $\mathbf{p}^{\ast}$ at time $1$. Each curve $\gamma\in\Gamma_{\delta}(t)$
is composed of a finite number of smooth curves, indexed by $k$. The $k$-th smooth curve $\gamma_{k}$
is defined on some time interval $[s_{k},s_{k+1}]\subset[0,1]$.
Then we have

$$ $d(t)=\inf_{\gamma=(\gamma_{k})_{k}\in\Gamma_{\delta}(t)}\sum_{k}\int_{s_{k}}^{ s_{k+1}}\|\gamma_{k}^{\prime}(s)\|_{\gamma_{k}(s)}ds,$ $$

where $\|\cdot\|_{\mathbf{p}}$ is the Riemmanian metric at $\mathbf{p}\in\mathcal{M}$.
As the closure of the open set $S_{\delta}$ is compact
and $\bar{S}_{\delta}$ is convex,
there are two positive constants $0<A_{\epsilon,\delta}\leq B_{\epsilon,\delta}$ such that

$$ $A_{\epsilon,\delta}\|\varphi(\mathbf{p}(t))-\varphi(\mathbf{p}^{\ast})\|_{ \epsilon}\leq d(t)\leq B_{\epsilon,\delta}\|\varphi(\mathbf{p}(t))-\varphi( \mathbf{p}^{\ast})\|_{\epsilon}.$ (38) $$

A relevant proof on how to obtain $A_{\epsilon,\delta}$ and $B_{\epsilon,\delta}$
is given in Hu ()[Lemma 3.3].

From (F), we have

$$ $\|\varphi(\mathbf{p}(t))-\varphi(\mathbf{p}^{\ast})\|_{\epsilon}\leq(\rho+ \epsilon)^{t+1}\|\varphi(\mathbf{p}(0))-\varphi(\mathbf{p}^{\ast})\|_{\epsilon}.$ (39) $$

Combining (38) and (39),
we obtain the constant $C=B_{\epsilon,\delta}/A_{\epsilon,\delta}$ such that

$$ $d(t)\leq C(\rho+\epsilon)^{t+1}d(0).$ $$

### F.1 Spectral radius of Jacobian matrix in the local cooridate

We next compute the matrix $\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$,
and then relate it to $\mathbf{T}^{\prime}(x^{\ast},y^{\ast})$.
This computation also tells us that
the tangent map of $\mathbf{T}$ at $(x^{\ast},y^{\ast})$
equals to $\mathbf{T}^{\prime}(x^{\ast},y^{\ast})$.

First of all, we specify the induced dynamics $\bar{\mathbf{T}}$ of
$\mathbf{T}$
by using the local coordinate representation of the retractions $\mathcal{R}_{1}$
and $\mathcal{R}_{2}$ near the fixed point $(x^{\ast},y^{\ast})$.
Let $(O_{1}\times O_{2},\varphi_{1}\times\varphi_{2})$ be the local chart.
Since a retraction (Absil et al., , Section 4.1)
is a smooth function from
the tangent bundle of a manifold to the manifold itself,
we can identify an open subset $B_{1}$ of the tangent bundle of
$\mathcal{M}_{1}$ (similarly on $B_{2}$ for $\mathcal{M}_{2}$) such that $x\in O_{1}$ and $\mathcal{R}_{1,x}(\xi_{1})\in O_{1}$
if $(x,\xi_{1})\in B_{1}$
This set $B_{1}$ can then be mapped to an open set $\bar{B}_{1}\subset\mathbb{R}^{d_{1}}\times\mathbb{R}^{d_{1}}$,
using the local chart of the tangent bundle ((Absil et al., , Section 3.5.3)).
On this subset $\bar{B}_{1}$,
we can define the local representation of $\mathcal{R}_{1}$ as follows

$$ $\bar{\mathcal{R}}_{1}:\bar{B}_{1}\to\bar{O}_{1},\quad\bar{\mathcal{R}}_{1}(u_{ 1},\bar{\xi}_{1})=\varphi_{1}\circ\mathcal{R}_{1,\varphi^{-1}_{1}(u_{1})}(D \varphi^{-1}_{1}(u_{1})[\bar{\xi}_{1}]),$ $$

Similarly for $\mathcal{R}_{2}$, we can define $\bar{\mathcal{R}}_{2}:\bar{B}_{2}\to\bar{O}_{2}$
on an open set $\bar{B}_{2}\subset\mathbb{R}^{d_{2}}\times\mathbb{R}^{d_{2}}$.
From the above definition of $\bar{\mathcal{R}}_{1}$ and $\bar{\mathcal{R}}_{2}$,
as well as the local coordinate of $\xi_{1}(x,y)$ and
$\xi_{2}(x,y)$ (as in (17),(18)),
we obtain the induced dynamics

$$ $\bar{\mathbf{T}}(u_{1},u_{2})=\left(\begin{array}[]{l}\bar{\mathcal{R}}_{1}(u_ {1},\bar{\xi}_{1}(u_{1},u_{2}))\\ \bar{\mathcal{R}}_{2}(u_{2},\bar{\xi}_{2}(u_{1},u_{2}))\end{array}\right),$ (40) $$

which is well defined on the non-empty open set
$\{(u_{1},u_{2})\in\bar{O}_{1}\times\bar{O}_{2}|(u_{1},\bar{\xi}_{1}(u_{1},u_{2
}))\in\bar{B}_{1},(u_{2},\bar{\xi}_{2}(u_{1},u_{2}))\in\bar{B}_{2}\}$.

From (40), we can compute the Jacobian matrix
$\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$ using the chain rule.
As $\bar{\xi}_{1}(u_{1}^{\ast},u_{2}^{\ast})=0$ and
$\bar{\xi}_{2}(u_{1}^{\ast},u_{2}^{\ast})=0$,
we have

$$ $\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})=\left(\begin{array}[]{ll} I_{d_{1}}+\frac{\partial\bar{\xi}_{1}}{\partial u_{1}}(u_{1}^{\ast},u_{2}^{ \ast})&\frac{\partial\bar{\xi}_{1}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast}) \\ \frac{\partial\bar{\xi}_{2}}{\partial u_{1}}(u_{1}^{\ast},u_{2}^{\ast})&I_{d_{ 2}}+\frac{\partial\bar{\xi}_{2}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast}) \end{array}\right),$ (41) $$

since the retraction by definition satisfies $\frac{\partial\bar{\mathcal{R}}_{1}}{\partial u_{1}}(u_{1}^{\ast},0)=I_{d_{1}}
,\frac{\partial\bar{\mathcal{R}}_{1}}{\partial\bar{\xi}_{1}}(u_{1}^{\ast},0)=I
_{d_{1}}$ (similarly for $\bar{\mathcal{R}}_{2}$).

We next verify that eigenvalues of the matrix in (41) are
the same as $\mathbf{T}^{\prime}(x^{\ast},y^{\ast})$.
Let $\delta^{\ast}=\sum_{i}\bar{\delta}_{i}^{\ast}E_{1,i}(x^{\ast})\in T_{x^{\ast}}
\mathcal{M}_{1}$,
$\eta^{\ast}=\sum_{j}\bar{\eta}_{j}^{\ast}E_{2,j}(y^{\ast})\in T_{y^{\ast}}
\mathcal{M}_{2}$,
then following the same argument as in (24) and (29), we have

$$ $\displaystyle\nabla_{x}\xi_{1}(x^{\ast},y^{\ast})[\delta^{\ast}]$ $\displaystyle=\sum_{i=1}^{d_{1}}\left(\frac{\partial\bar{\xi}_{1,i}}{\partial u _{1}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\delta}^{\ast}\right)E_{1,i}(x^{\ast}),$ $\displaystyle D_{y}\xi_{1}(x^{\ast},y^{\ast})[\eta^{\ast}]$ $\displaystyle=\sum_{i=1}^{d_{1}}\left(\frac{\partial\bar{\xi}_{1,i}}{\partial u _{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}\right)E_{1,i}(x^{\ast}),$ $\displaystyle D_{x}\xi_{2}(x^{\ast},y^{\ast})[\delta^{\ast}]$ $\displaystyle=\sum_{j=1}^{d_{2}}\left(\frac{\partial\bar{\xi}_{2,j}}{\partial u _{1}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\delta}^{\ast}\right)E_{2,j}(y^{\ast}),$ $\displaystyle\nabla_{y}\xi_{2}(x^{\ast},y^{\ast})[\eta^{\ast}]$ $\displaystyle=\sum_{j=1}^{d_{2}}\left(\frac{\partial\bar{\xi}_{2,j}}{\partial u _{2}}(u_{1}^{\ast},u_{2}^{\ast})\bar{\eta}^{\ast}\right)E_{2,j}(y^{\ast}).$ $$

From the definition of eigenvalue and eigenvector pairs,
these four equations indicate that their eigenvalues are indeed the same.
They also indicate that the tangent map of $\mathbf{T}$ at $(x^{\ast},y^{\ast})$
equals to (9),
which means that the tangent map of $\mathbf{T}$ at the fixed point
does not depend on the choice of the retraction.

## Appendix G Proof of Theorem 3.2

We first rewrite $\mathbf{T}^{\prime}(x^{\ast},y^{\ast})=I+\gamma\mathbf{M}_{g}$
in a local coordinate chart with
$\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})=I_{d_{1}+d_{2}}+\gamma M_
{g}^{\prime}$.

From the definition of $\mathbf{M}_{g}$ in (11)
and the connection between $\mathbf{T}^{\prime}$ and $\bar{\mathbf{T}}^{\prime}$
in (41),(31)-(33), we have

$$ $M_{g}^{\prime}=\left(\begin{array}[]{ll}-\bar{g}_{1}(u_{1}^{\ast})^{-1}\cdot C &-\bar{g}_{1}(u_{1}^{\ast})^{-1}\cdot B\\ \tau\bar{g}_{2}(u_{2}^{\ast})^{-1}\cdot B^{\intercal}&-\tau\bar{g}_{2}(u_{2}^{ \ast})^{-1}\cdot A\end{array}\right),$ $$

where the matrices $A$,$B$,$C$ are defined in (26).
Furthermore, $\mathbf{M}_{g}$ and $M_{g}^{\prime}$ have the same eigenvalues.

We next show that the real part of each eigenvalue of $M_{g}^{\prime}$ is strictly smaller than zero.
From this, $\gamma^{\bullet}(\mathbf{M}_{g})=\gamma^{\bullet}(M_{g}^{\prime})>0$.
It is easy to check that if $0<\gamma<\gamma^{\bullet}(M_{g}^{\prime})$,
the spectral radius of $\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$ is strictly smaller than one.

To analyze $M_{g}^{\prime}$, we could not apply (Li et al., , Lemma 5.2)
directly since the local metric $\bar{g}_{1}(u_{1}^{\ast})$
and $\bar{g}_{2}(u_{2}^{\ast})$ are not identity matrices. We consider a matrix
which is similar to $M_{g}^{\prime}$ so that they have the same eigenvalues,

$$ $M_{g}=\left(\begin{array}[]{ll}-\bar{\mathbf{C}}&-\bar{\mathbf{B}}\\ \tau\bar{\mathbf{B}}^{\intercal}&-\tau\bar{\mathbf{A}}\end{array}\right),$ (42) $$

where

- •
$\bar{\mathbf{C}}=\bar{g}_{1}(u_{1}^{\ast})^{-1/2}\cdot C\cdot\bar{g}_{1}(u_{1}
^{\ast})^{-1/2}$
- •
$\bar{\mathbf{B}}=\bar{g}_{1}(u_{1}^{\ast})^{-1/2}\cdot B\cdot\bar{g}_{2}(u_{2}
^{\ast})^{-1/2}$
- •
$\bar{\mathbf{A}}=\bar{g}_{2}(u_{2}^{\ast})^{-1/2}\cdot A\cdot\bar{g}_{2}(u_{2}
^{\ast})^{-1/2}$

Under the assumptions of Theorem 3.2,
we check that the following conditions hold:

- •
$\bar{\mathbf{A}}$ is p.d because $\mathbf{A}$ is p.d.
- •
$\bar{\mathbf{C}}+\bar{\mathbf{B}}\bar{\mathbf{A}}^{-1}\bar{\mathbf{B}}^{\intercal}$ is p.d because $\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B}^{\intercal}$ is p.d.
- •
$\tau>\|\bar{\mathbf{C}}\|/\lambda_{min}(\bar{\mathbf{A}})$.

Indeed, we can relate these conditions to the following intrinsic quantities on the manifold $\mathcal{M}_{1}\times\mathcal{M}_{2}$, by following the proof of Theorem LABEL:thm:locdse (see Appendix A):

- •
Let $\mathbf{A}=-\mbox{Hess}_{y}f(x^{\ast},y^{\ast})$, then $\mathbf{A}$ and $\bar{\mathbf{A}}$ have the same eigenvalues.
- •
Let $\mathbf{B}=\mbox{grad}^{2}_{yx}f(x^{\ast},y^{\ast})$, $\mathbf{B}^{\intercal}=\mbox{grad}^{2}_{xy}f(x^{\ast},y^{\ast})$, then $\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B}^{\intercal}$ and $\bar{\mathbf{C}}+\bar{\mathbf{B}}\bar{\mathbf{A}}^{-1}\bar{\mathbf{B}}^{\intercal}$ have the same eigenvalues.
- •
Let $\mathbf{C}=\mbox{Hess}_{x}f(x^{\ast},y^{\ast})$, then $\mathbf{C}$ and $\bar{\mathbf{C}}$ have the same eigenvalues. As $\mathbf{C}$ and $\bar{\mathbf{C}}$ are symmetric, $\|\mathbf{C}\|=\max_{k}|\lambda_{k}(\mathbf{C})|=\|\bar{\mathbf{C}}\|$. Therefore,
the condition $\tau>\|\mathbf{C}\|/\lambda_{min}(\mathbf{A})$ is equivalent to $\tau>\|\bar{\mathbf{C}}\|/\lambda_{min}(\bar{\mathbf{A}})$.

### G.1 Spectral analysis of M g subscript 𝑀 𝑔 M_{g} italic_M start_POSTSUBSCRIPT italic_g end_POSTSUBSCRIPT

To analyze the eigenvalues of $M_{g}$, we
use the next proposition
which is adapted from (Li et al., , Lemma 5.2)
with a refined range of $\tau$.
For two real symmetric matrices $A$ and $B$,
we write $A\geq B$ if $A-B$ is semi-p.d.
We write $A>B$ if $A-B$ is p.d.

Let us first introduce a working assumption which will be needed several times.

###### Assumption G.1 .

Let $A\in\mathbb{R}^{m\times m}$, $B\in\mathbb{R}^{n\times m}$, $C\in\mathbb{R}^{n\times n}$
such that $A$ and $C+BA^{-1}B^{\intercal}$ are positive definite.

###### Proposition G.1 .

Under Assumption G.1 and $\tau>\|C\|/\lambda_{min}(A)$,
any eigenvalue $\lambda=\lambda_{0}+i\lambda_{1}$ (with $\lambda_{0}\in\mathbb{R},\lambda_{1}\in\mathbb{R}$) of the following matrix,

$$ $M=\left(\begin{array}[]{ll}-C&-B\\ \tau B^{\intercal}&-\tau A\end{array}\right)$ $$

satisfies $\lambda_{0}<0$.

###### Proof.

We show that $\lambda_{0}\geq 0$ will lead to a contradiction,
by following the proof of (Li et al., , Lemma 5.2).
Assume that $\lambda$ is an eigenvalue of $M$, i.e. $\mbox{det}(\lambda I-M)=0$.
To compute $\mbox{det}(\lambda I-M)$, note that
$\lambda I+\tau A$ is invertible since
$\lambda_{0}I+\tau A$ is positive definite ($\tau A$ is positive definite and $\lambda_{0}\geq 0$).
By the Schur complement,
we have

$$ $\displaystyle\mbox{det}(\lambda I-M)$ $\displaystyle=\mbox{det}\left(\begin{array}[]{cc}C+\lambda I&B\\ -\tau B^{\intercal}&\tau A+\lambda I\end{array}\right)$ $\displaystyle=\mbox{det}(\lambda I+\tau A)\mbox{det}(H(\lambda))$ $$

where $H(\lambda)=\lambda I+C+B(\lambda/\tau I+A)^{-1}B^{\intercal}$.
As $\lambda I+\tau A$ is invertible, $\mbox{det}(\lambda I-M)=0$ implies that
$\mbox{det}(H(\lambda))=0$.

Let the spectral decomposition of $A$ be $U\Lambda_{A}U^{\intercal}$, with orthogonal $U\in\mathbb{R}^{m\times m}$ and diagonal $\Lambda_{A}=\mbox{diag}(\lambda_{1}(A),\cdots,\lambda_{m}(A))\in\mathbb{R}^{m
\times m}$.
Then

$$ $H(\lambda)=\lambda I+C+\tilde{B}D\tilde{B}^{\intercal},$ $$

with $\tilde{B}=BU$ and $D=\mbox{diag}(d_{1},\cdots,d_{m})$ where

$$ $d_{k}=\frac{1}{\lambda/\tau+\lambda_{k}(A)}=\frac{\lambda_{0}/\tau+\lambda_{k} (A)-i\lambda_{1}/\tau}{(\lambda_{0}/\tau+\lambda_{k}(A))^{2}+(\lambda_{1}/\tau )^{2}},\quad 1\leq k\leq m.$ $$

It follows that if $\lambda_{0}>\|C\|$, the real-part of $H(\lambda)$ is

$$ $\mbox{Re}(H(\lambda))=\lambda_{0}I+C+\tilde{B}\mbox{Re}(D)\tilde{B}^{\intercal }\quad\mbox{p.d}.$ (43) $$

This is contradictory to the fact that $\mbox{det}(H(\lambda))=0$ (Li et al., , Corollary 10.2).

On the other hand, if $0\leq\lambda_{0}\leq\|C\|$ and $\lambda_{1}\neq 0$, we consider for $\beta\in\mathbb{R}$,

$$ $\displaystyle\mbox{Re}(H(\lambda))+\frac{\tau\beta}{\lambda_{1}}\mbox{Im}(H( \lambda))$ $\displaystyle=\lambda_{0}I+C+\tilde{B}\mbox{Re}(D)\tilde{B}^{\intercal}+\frac{ \tau\beta}{\lambda_{1}}(\lambda_{1}I+\tilde{B}\mbox{Im}(D)\tilde{B}^{\intercal})$ $\displaystyle=(\lambda_{0}+\tau\beta)I+C+\tilde{B}F\tilde{B}^{\intercal},$ (44) $$

where
$F=\mbox{diag}(f_{1},\cdots,f_{m})$ with

$$ $\displaystyle f_{k}=\frac{\lambda_{0}/\tau+\lambda_{k}(A)-\beta}{(\lambda_{0}/ \tau+\lambda_{k}(A))^{2}+(\lambda_{1}/\tau)^{2}},\quad 1\leq k\leq m.$ (45) $$

Take $\beta=\lambda_{min}(A)$, then $f_{k}\geq 0$ for each $k\leq m$.
The condition $\tau\lambda_{min}(A)>\|C\|$ implies that $(\lambda_{0}+\tau\beta)I+C$ is p.d.
and together with (45),
we have (44) is p.d., so it is contradictory to
the fact that $\mbox{det}(H(\lambda))=0$ (Li et al., , Lemma 10.1).

Lastly, if $0\leq\lambda_{0}\leq\|C\|$ and $\lambda_{1}=0$, we have from (43),

$$ $H(\lambda)=\lambda_{0}I+C+\tilde{B}D\tilde{B}^{\intercal}$ (46) $$

with $d_{k}=\frac{1}{\lambda_{0}/\tau+\lambda_{k}(A)}=\frac{1}{\lambda_{k}(A)}-\frac
{\lambda_{0}/\tau}{(\lambda_{0}/\tau+\lambda_{k}(A))\lambda_{k}(A)}\geq\frac{1
}{\lambda_{k}(A)}-\frac{\lambda_{0}/\tau}{\lambda_{min}(A)\lambda_{k}(A)}$.
It follows that

$$ $\displaystyle H(\lambda)$ $\displaystyle\geq\lambda_{0}I+C+\bigg{(}1-\frac{\lambda_{0}/\tau}{\lambda_{min }(A)}\bigg{)}BA^{-1}B^{\intercal}$ $\displaystyle=\bigg{(}1-\frac{\lambda_{0}/\tau}{\lambda_{min}(A)}\bigg{)}(C+BA ^{-1}B^{\intercal})+\frac{\lambda_{0}/\tau}{\lambda_{min}(A)}C+\lambda_{0}I.$ $$

As $\tau>\|C\|/\lambda_{min}(A)$, $0\leq\lambda_{0}\leq\|C\|$ we have that $0\leq\frac{\lambda_{0}/\tau}{\lambda_{min}(A)}<1$ and that $I+\frac{1}{\tau\lambda_{min}(A)}C$ is p.d,
i.e. we find that again $H(\lambda)$ is p.d which is contradictory. In conclusion, $\lambda_{0}<0$.
∎

###### Assumption G.1 .

###### Proposition G.1 .

###### Proof.

### G.2 Local convergence rate of τ 𝜏 \tau italic_τ -GDA

We use the next result obtained in Li et al. (, Lemma 5.3)
to control the local convergence rate of $\tau$-GDA.
It controls the spectral radius of
the matrix $I+\gamma M_{g}$
on a specific range of $\tau$ and $\gamma$.

Recall that $L_{g}=\max(\|\mathbf{A}\|,\|\mathbf{B}\|,\|\mathbf{C}\|)$
and $\mu_{g}=\min(L_{g},\lambda_{min}(\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B
}^{\intercal}))$.

###### Proposition G.2 .

Assume $(x^{\ast},y^{\ast})$ is a DSE of $f\in C^{2}$.
If $\tau\geq\frac{2L_{g}}{\lambda_{min}(\mathbf{A})}$
and $\gamma=\frac{1}{4\tau L_{g}}$, we have
$\rho(I+\gamma M_{g})\leq 1-\frac{\mu_{g}}{16\tau L_{g}}$.

###### Proof.

Let $L=\max(\|\bar{\mathbf{A}}\|,\|\bar{\mathbf{B}}\|,\|\bar{\mathbf{C}}\|)$
and $\mu_{x}=\min(L,\lambda_{min}(\bar{\mathbf{C}}+\bar{\mathbf{B}}\bar{\mathbf{A}}
^{-1}\bar{\mathbf{B}}^{\intercal}))$.
We verify that $\lambda_{min}(\mathbf{A})=\lambda_{\min}(\bar{\mathbf{A}})$, $L=L_{g}$ and $\mu_{g}=\mu_{x}$.
From the proof of Li et al. (, Lemma 5.3),
we obtain directly the upper bound of the spectral radius of $I+\gamma M_{g}$,
related to $L$ and $\mu_{x}$.
∎

###### Proposition G.2 .

###### Proof.

## Appendix H Proof of Theorem 3.3

Based on the proof of Jin et al. (, Proposition 26),
we consider a
matrix $\tilde{M}_{g}$ which is similar to $M_{g}$ (defined in (42)).
We shall verify that for any $\tau>0$,
the eigenvalues of $\tilde{M}_{g}$ are strictly negative, therefore
one can compute $\gamma^{\bullet}(\mathbf{M}_{g})=\gamma^{\bullet}(M_{g})=\gamma^{\bullet}(
\tilde{M}_{g})$
to set the range for $\gamma$.

$$ $\displaystyle\tilde{M}_{g}$ $\displaystyle=\left(\begin{array}[]{ll}0&\sqrt{\frac{1}{\tau}}I\\ I&0\end{array}\right)\left(\begin{array}[]{ll}-\bar{\mathbf{C}}&-\bar{\mathbf{ B}}\\ \tau\bar{\mathbf{B}}^{\intercal}&-\tau\bar{\mathbf{A}}\end{array}\right)\left( \begin{array}[]{ll}0&\sqrt{\frac{1}{\tau}}I\\ I&0\end{array}\right)^{-1}$ $\displaystyle=\left(\begin{array}[]{ll}\sqrt{\tau}\bar{\mathbf{B}}^{\intercal} &-\sqrt{\tau}\bar{\mathbf{A}}\\ -\bar{\mathbf{C}}&-\bar{\mathbf{B}}\end{array}\right)\left(\begin{array}[]{ll} 0&I\\ \sqrt{\tau}I&0\end{array}\right)=\left(\begin{array}[]{ll}-\tau\bar{\mathbf{A} }&\sqrt{\tau}\bar{\mathbf{B}}^{\intercal}\\ -\sqrt{\tau}\bar{\mathbf{B}}&-\bar{\mathbf{C}}\end{array}\right).$ $$

From the Ky Fan inequality,
the real part of each eigenvalue $\lambda$ of $\tilde{M}_{g}$
is upper bounded by the maximal eigenvalue of
$(\tilde{M}_{g}+\tilde{M}_{g}^{\intercal})/2$.
We verify that indeed $\lambda_{max}((\tilde{M}_{g}+\tilde{M}_{g}^{\intercal})/2)<0$
because $\bar{\mathbf{A}}$ and $\bar{\mathbf{C}}$ are p.d., and
$\tau>0$.
Therefore, the real part of $\lambda$ is strictly smaller than 0 0.

## Appendix I Derivation of deterministic τ 𝜏 \tau italic_τ -SGA algorithm

The $\tau$-SGA algorithm modifies the vector field
$\xi(x,y)=(-\delta(x,y),\tau\eta(x,y))$
of $\tau$-GDA using the anti-symmetric part of
the Jacobian matrix of $\xi(x,y)$.
As $\delta(x,y)=\mbox{grad}_{x}f(x,y)$ and
$\eta(x,y)=\mbox{grad}_{y}f(x,y)$,
the Jacobian matrix is the following linear transform
on $T_{x}\mathcal{M}_{1}\times T_{y}\mathcal{M}_{2}$
(which is a natural extension from the Euclidean case),

$$ $J(x,y)=\left(\begin{array}[]{ll}-\tilde{\mathbf{C}}(x,y)&-\tilde{\mathbf{B}}(x ,y)\\ \tau\tilde{\mathbf{B}}^{\intercal}(x,y)&-\tau\tilde{\mathbf{A}}(x,y)\end{array }\right)=\left(\begin{array}[]{ll}-\mbox{Hess}_{x}(x,y)&-\mbox{grad}^{2}_{yx}f (x,y)\\ \tau\mbox{grad}^{2}_{xy}f(x,y)&\tau\mbox{Hess}_{y}(x,y)\end{array}\right).$ (47) $$

Here we introduce the notation $\tilde{\mathbf{A}},\tilde{\mathbf{B}},\tilde{\mathbf{C}}$ in (47) to simplify the equation.
The $\tau$-SGA update rule is obtained from
(^5^55 Note that in the original SGA rule (Letcher et al., , Proposition 5) the transpose of $\frac{J(x,y)-J^{\intercal}(x,y)}{2}$ is considered since their definition of $J$ has a sign difference compared to the $J$ in (47).)

$$ $\displaystyle\xi(x,y)+\mu\left(\frac{J(x,y)-J^{\intercal}(x,y)}{2}\right)\xi(x ,y)$ $\displaystyle=\xi(x,y)+\mu\frac{\tau+1}{2}\left(\begin{array}[]{ll}0&-\tilde{ \mathbf{B}}(x,y)\\ \tilde{\mathbf{B}}^{\intercal}(x,y)&0\end{array}\right)\xi(x,y)$ $\displaystyle=\left[\left(\begin{array}[]{l}-\delta\\ \tau\eta\end{array}\right)+\mu\frac{\tau+1}{2}\left(\begin{array}[]{l}-\tau \tilde{\mathbf{B}}[\eta]\\ -\tilde{\mathbf{B}}^{\intercal}[\delta]\end{array}\right)\right](x,y).$ $$

We next show that the correction term, which is proportional to $(\tau\tilde{\mathbf{B}}[\eta],\tilde{\mathbf{B}}^{\intercal}[\delta])$,
is orthogonal to the $\tau$-GDA direction $(-\delta,\tau\eta)$
at each $(x,y)$,
under the Riemannian metric on the tangent space $T_{x}\mathcal{M}_{1}\times T_{y}\mathcal{M}_{2}$.

### I.1 Proof of Proposition 3.1

We follow the proof in Appendix A,
by using
a local coordinate chart $(O_{1}\times O_{2},\varphi_{1}\times\varphi_{2})$ around the point $(x,y)$ rather than $(x^{\ast},y^{\ast})$.
As in (17) and (18), this
coordinate chart maps $\delta(x,y)$ and $\eta(x,y)$ to
their local coordinates $\bar{\delta}(u_{1},u_{2})$
and $\bar{\eta}(u_{1},u_{2})$. It also induces a canonical basis
$\{E_{1,i}(x)\}_{i\leq d_{1}}$ on $T_{x}\mathcal{M}_{1}$ and $\{E_{2,j}(y)\}_{j\leq d_{2}}$ on $T_{y}\mathcal{M}_{2}$ .
By the definition of cross-gradients,

$$ $\displaystyle\tilde{\mathbf{B}}(x,y)[\eta(x,y)]$ $\displaystyle=D_{y}\mbox{grad}_{x}f(x,y)[\eta(x,y)]$ $\displaystyle=D_{y}\delta(x,y)[\eta(x,y)]$ $\displaystyle=\sum_{i}\left(\frac{\partial\bar{\delta}_{i}}{\partial u_{2}}(u_ {1},u_{2})\bar{\eta}(u_{1},u_{2})\right)E_{1,i}(x).$ $$

Using the Riemannian metric $\bar{g}_{1}$ represented in the local coordinate
as in (19), it turns out that

$$ $\langle\tilde{\mathbf{B}}(x,y)[\eta(x,y)],\delta(x,y)\rangle_{x}=\bar{\delta}( u_{1},u_{2})^{\intercal}\bar{g}_{1}(u_{1})\left(\frac{\partial\bar{\delta}}{ \partial u_{2}}(u_{1},u_{2})\bar{\eta}(u_{1},u_{2})\right).$ (48) $$

Similarly we have

$$ $\langle\tilde{\mathbf{B}}^{\intercal}(x,y)[\delta(x,y)],\eta(x,y)\rangle_{y}= \bar{\eta}(u_{1},u_{2})^{\intercal}\bar{g}_{2}(u_{2})\left(\frac{\partial\bar{ \eta}}{\partial u_{1}}(u_{1},u_{2})\bar{\delta}(u_{1},u_{2})\right).$ (49) $$

We conclude that
(48) equals to (49) because
as in (32),(33)

$$ $\bar{g}_{1}(u_{1})\frac{\partial\bar{\delta}}{\partial u_{2}}(u_{1},u_{2})= \left(\bar{g}_{2}(u_{2})\frac{\partial\bar{\eta}}{\partial u_{1}}(u_{1},u_{2}) \right)^{\intercal}.$ $$

## Appendix J Proof of Theorem 3.4

Since $f$ is twice continuously differentiable,
we apply Theorem 3.1
to analyze the induced dynamics $\bar{\mathbf{T}}$ (where $\mathbf{T}$
is the update rule of Asymptotic $\tau$-SGA).
Following (40),
we obtain

$$ $\bar{\mathbf{T}}(u_{1},u_{2})=\left(\begin{array}[]{l}\bar{\mathcal{R}}_{1} \left(u_{1},-\gamma\left[\bar{\delta}+\mu\frac{(\tau+1)\tau}{2}\frac{\partial \bar{\delta}}{\partial u_{2}}\bar{\eta}\right](u_{1},u_{2})\right)\\ \bar{\mathcal{R}}_{2}\left(u_{2},\gamma\tau\bar{\eta}(u_{1},u_{2})\right)\end{ array}\right).$ (50) $$

From (50), we compute the Jacobian matrix
$\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})=I+\gamma M_{s}^{\prime}$,
where

$$ $M_{s}^{\prime}=\left(\begin{array}[]{ll}-\frac{\partial\bar{\delta}}{\partial u _{1}}(u_{1}^{\ast},u_{2}^{\ast})&-\frac{\partial\bar{\delta}}{\partial u_{2}}( u_{1}^{\ast},u_{2}^{\ast})\\ \tau\frac{\partial\bar{\eta}}{\partial u_{1}}(u_{1}^{\ast},u_{2}^{\ast})&\tau \frac{\partial\bar{\eta}}{\partial u_{2}}(u_{1}^{\ast},u_{2}^{\ast})\end{array }\right)-\mu\frac{(\tau+1)\tau}{2}\left(\begin{array}[]{ll}[\frac{\partial\bar {\delta}}{\partial u_{2}}\frac{\partial\bar{\eta}}{\partial u_{1}}](u_{1}^{ \ast},u_{2}^{\ast})&[\frac{\partial\bar{\delta}}{\partial u_{2}}\frac{\partial \bar{\eta}}{\partial u_{2}}](u_{1}^{\ast},u_{2}^{\ast})\\ 0&0\end{array}\right).$ $$

This can be reduced to the analysis of the eigenvalues of a similar matrix $M_{s}$
as in (42),

$$ $M_{s}=\left(\begin{array}[]{ll}-\bar{\mathbf{C}}&-\bar{\mathbf{B}}\\ \tau\bar{\mathbf{B}}^{\intercal}&-\tau\bar{\mathbf{A}}\end{array}\right)+\mu \frac{(\tau+1)\tau}{2}\left(\begin{array}[]{ll}-\bar{\mathbf{B}}\bar{\mathbf{B }}^{\intercal}&\bar{\mathbf{B}}\bar{\mathbf{A}}\\ 0&0\end{array}\right).$ $$

As in the proof of Theorem 3.2,
we verify that $M_{s}$ and the $\mathbf{M}_{s}$ in (14)
have the same eigenvalues.
From the assumption of Theorem 3.4, we verify that
$\tau>\min(\|\bar{\mathbf{C}}\|,\|\bar{\mathbf{C}}+\theta\bar{\mathbf{B}}\bar{
\mathbf{B}}^{\intercal}\|)/\lambda_{min}(\bar{\mathbf{A}})$,
and $0\leq\theta\leq 1/\lambda_{max}(\bar{\mathbf{A}})$.
We can therefore apply
Proposition J.1 (see next)
to conclude that the real-part of each eigenvalue of $\mathbf{M}_{s}$ is strictly negative.
Therefore for $0<\gamma<\gamma^{\bullet}(\mathbf{M}_{s})$,
Asymptotic $\tau$-SGA is locally convergent to DSE
with rate $\rho(I+\gamma\mathbf{M}_{s})$.

Before we proceed to analyze $M_{s}$, we remark that
the non-asymptotic analysis of $\tau$-SGA remains an interesting open question.
Indeed, if we want to analyze $\tau$-SGA beyond
the asymptotic regime (e.g. $\tau$ is small),
one needs to consider this induced dynamics
$\bar{\mathcal{R}}_{2}\left(u_{2},\gamma\left[\tau\bar{\eta}-\mu\frac{\tau+1}{2
}\frac{\partial\bar{\eta}}{\partial u_{1}}\bar{\delta}\right](u_{1},u_{2})\right)$ for the variable $u_{2}$.
The analysis of the spectral radius of the Jacobian matrix
$\bar{\mathbf{T}}^{\prime}(u_{1}^{\ast},u_{2}^{\ast})$ is harder as one
could not easily adapt the proof of Proposition J.1 to this case.

### J.1 Spectral analysis of M s subscript 𝑀 𝑠 M_{s} italic_M start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT

The following proposition is adapted from Proposition G.1
to analyze the eigenvalues of $M_{s}$.

###### Proposition J.1 .

Under Assumption G.1,
$\tau>\frac{\min(\|C\|,\|C+\theta BB^{\intercal}\|)}{\lambda_{min}(A)}$,
and $\mu=\theta\frac{2}{\tau(\tau+1)}$ with $0\leq\theta\leq 1/\lambda_{max}(A)$,
any eigenvalue $\lambda=\lambda_{0}+i\lambda_{1}$ (where $\lambda_{0}\in\mathbb{R},\lambda_{1}\in\mathbb{R}$) of the following matrix

$$ $M=\left(\begin{array}[]{ll}-C&-B\\ \tau B^{\intercal}&-\tau A\end{array}\right)+\mu\frac{(\tau+1)\tau}{2}\left( \begin{array}[]{ll}-BB^{\intercal}&BA^{\intercal}\\ 0&0\end{array}\right)$ $$

satisfies $\lambda_{0}<0$.

###### Proof.

We show that $\lambda_{0}\geq 0$ will lead to a contradiction.
Assume that $\lambda$ is an eigenvalue of $M$, i.e. $\mbox{det}(\lambda I-M)=0$.
By following the proof of Proposition G.1,
we have

$$ $\displaystyle\mbox{det}(\lambda I-M)$ $\displaystyle=\mbox{det}\left(\begin{array}[]{cc}C+\theta BB^{\intercal}+ \lambda I&B-\theta BA\\ -\tau B^{\intercal}&\tau A+\lambda I\end{array}\right)$ (53) $\displaystyle=\mbox{det}(\lambda I+\tau A)\mbox{det}(H(\lambda))$ (54) $$

where $H(\lambda)=\lambda I+C+\theta BB^{\intercal}+(B-\theta BA)(\lambda/\tau I+A)^{
-1}B^{\intercal}$.
As $\lambda I+\tau A$ is invertible, $\mbox{det}(\lambda I-M)=0$ implies that
$\mbox{det}(H(\lambda))=0$.

Let the spectral decomposition of $A$ be $U\Lambda_{A}U^{\intercal}$, with orthogonal $U\in\mathbb{R}^{m\times m}$ and diagonal $\Lambda_{A}=\mbox{diag}(\lambda_{1}(A),\cdots,\lambda_{m}(A))\in\mathbb{R}^{m
\times m}$.
Then

$$ $H(\lambda)=\lambda I+C+\theta BB^{\intercal}+\tilde{B}D\tilde{B}^{\intercal},$ (55) $$

with $\tilde{B}=BU$ and $D=\mbox{diag}(d_{1},\cdots,d_{m})$ where

$$ $d_{k}=\frac{1-\theta\lambda_{k}(A)}{\lambda/\tau+\lambda_{k}(A)}=(1-\theta \lambda_{k}(A))\frac{\lambda_{0}/\tau+\lambda_{k}(A)-i\lambda_{1}/\tau}{( \lambda_{0}/\tau+\lambda_{k}(A))^{2}+(\lambda_{1}/\tau)^{2}},\quad 1\leq k\leq m.$ $$

It follows that if $\lambda_{0}>\min(\|C\|,\|C+\theta BB^{\intercal}\|)$,
the real-part of $H(\lambda)$ is

$$ $\mbox{Re}(H(\lambda))=\lambda_{0}I+C+\theta BB^{\intercal}+\tilde{B}\mbox{Re}( D)\tilde{B}^{\intercal}\geq\lambda_{0}I+C+\theta BB^{\intercal}\quad\mbox{p.d}$ (56) $$

This is contradictory to the fact that $\mbox{det}(H(\lambda))=0$ (Li et al., , Corollary 10.2).

On the other hand, if $0\leq\lambda_{0}\leq\min(\|C\|,\|C+\theta BB^{\intercal}\|)$ and $\lambda_{1}\neq 0$,
we consider for $\beta\in\mathbb{R}$,

$$ $\displaystyle\mbox{Re}(H(\lambda))+\frac{\tau\beta}{\lambda_{1}}\mbox{Im}(H( \lambda))$ $\displaystyle=\lambda_{0}I+C+\theta BB^{\intercal}+\tilde{B}\mbox{Re}(D)\tilde {B}^{\intercal}+\frac{\tau\beta}{\lambda_{1}}(\lambda_{1}I+\tilde{B}\mbox{Im}( D)\tilde{B}^{\intercal})$ $\displaystyle=(\lambda_{0}+\tau\beta)I+C+\theta BB^{\intercal}+\tilde{B}F \tilde{B}^{\intercal},$ (57) $$

where
$F=\mbox{diag}(f_{1},\cdots,f_{m})$ with

$$ $\displaystyle f_{k}=(1-\theta\lambda_{k}(A))\frac{\lambda_{0}/\tau+\lambda_{k} (A)-\beta}{(\lambda_{0}/\tau+\lambda_{k}(A))^{2}+(\lambda_{1}/\tau)^{2}},\quad 1 \leq k\leq m.$ (58) $$

Take $\beta=\lambda_{min}(A)$, then $f_{k}\geq 0$ for each $k\leq m$.
The condition $\tau\lambda_{min}(A)>\min(\|C\|,\|C+\theta BB^{\intercal}\|)$ implies
that $(\lambda_{0}+\tau\beta)I+C+\theta BB^{\intercal}$ is p.d.
and together with (58),
we get that (57) is p.d., so it is contradictory to
the fact that $\mbox{det}(H(\lambda))=0$ (Li et al., , Lemma 10.1).

Lastly, if $0\leq\lambda_{0}\leq\min(\|C\|,\|C+\theta BB^{\intercal}\|)$
and $\lambda_{1}=0$, we have from (56)

$$ $H(\lambda)=\lambda_{0}I+C+\theta BB^{\intercal}+\tilde{B}D\tilde{B}^{\intercal}$ $$

with $d_{k}=\frac{1-\theta\lambda_{k}(A)}{\lambda_{0}/\tau+\lambda_{k}(A)}=\frac{1-
\theta\lambda_{k}(A)}{\lambda_{k}(A)}-\frac{(1-\theta\lambda_{k}(A))\lambda_{0
}/\tau}{(\lambda_{0}/\tau+\lambda_{k}(A))\lambda_{k}(A)}\geq(1-\theta\lambda_{
k}(A))(\frac{1}{\lambda_{k}(A)}-\frac{\lambda_{0}/\tau}{\lambda_{min}(A)
\lambda_{k}(A)})$.
As $\tau>\min(\|C\|,\|C+\theta BB^{\intercal}\|)/\lambda_{min}(A)$
and $0\leq\lambda_{0}\leq\min(\|C\|,\|C+\theta BB^{\intercal}\|)$ we have $0\leq\frac{\lambda_{0}/\tau}{\lambda_{min}(A)}<1$ and it follows

$$ $\displaystyle H(\lambda)$ $\displaystyle\geq\lambda_{0}I+C+\theta BB^{\intercal}+\left(1-\frac{\lambda_{0 }/\tau}{\lambda_{min}(A)}\right)BA^{-1}B^{\intercal}-\theta(1-\frac{\lambda_{0 }/\tau}{\lambda_{min}(A)})\tilde{B}\tilde{B}^{\intercal}$ $\displaystyle=\left(1-\lambda_{0}c_{0}\right)(C+BA^{-1}B^{\intercal})+\lambda_ {0}(c_{0}C+c_{0}\theta BB^{\intercal}+I)$ $\displaystyle\geq\left(1-\lambda_{0}c_{0}\right)(C+BA^{-1}B^{\intercal})$ (59) $$

where $c_{0}=\frac{1}{\tau\lambda_{min}(A)}$. The last inequity (59) is
due to $\lambda_{0}\geq 0$ and the condition $1>c_{0}\min(\|C\|,\|C+\theta BB^{\intercal}\|)$.
They imply that $\lambda_{0}(c_{0}C+c_{0}\theta BB^{\intercal}+I)\geq 0$.
As $1-\lambda_{0}c_{0}\in(0,1]$, (59) implies
that $H(\lambda)$ is p.d which is contradictory. In conclusion, $\lambda_{0}<0$.
∎

###### Proposition J.1 .

###### Proof.

### J.2 Local convergence rate of Asymptotic τ 𝜏 \tau italic_τ -SGA

The local convergence analysis
in Proposition G.2
provides an upper bound of the rate $\rho(I+\gamma M_{g})$
of $\tau$-GDA.
This section extends this result
to Asymptotic $\tau$-SGA.
We aim to obtain an upper bound of the rate $\rho(I+\gamma M_{s})$
which is smaller than that of $\rho(I+\gamma M_{g})$.

The following lemma analyzes the eigenvalues of $M_{s}$ for Asymptotic $\tau$-SGA
when $\mu=\theta\frac{2}{\tau(\tau+1)}$.
It refines the eigenvalue bounds
of Li et al. (, Lemma 5.2)
on $M_{g}$ of $\tau$-GDA (when $\mu=0$).

###### Lemma J.1 .

Under Assumption G.1,
$\tau>0$, and $0\leq\theta\leq 1/\lambda_{max}(A)$,
any eigenvalue $\lambda=\lambda_{0}+i\lambda_{1}$ (where $\lambda_{0}\in\mathbb{R},\lambda_{1}\in\mathbb{R}$) of the following matrix

$$ $M=\left(\begin{array}[]{ll}-C&-B\\ \tau B^{\intercal}&-\tau A\end{array}\right)+\theta\left(\begin{array}[]{ll}- BB^{\intercal}&BA^{\intercal}\\ 0&0\end{array}\right)$ $$

satisfies

- a.
$|\lambda_{1}|\leq\sqrt{\tau}\sqrt{1-\theta\lambda_{min}(A)}\|B\|$.
- b.
If $\lambda_{1}\neq 0$, we have $\lambda_{0}\leq-\lambda_{+}$ with $\lambda_{+}=\frac{1}{2}(\lambda_{min}(A)\tau-\min(\|C\|,\|C+\theta BB^{
\intercal}\|))$.
- c.
Let $\tau=\frac{\min(\|C\|,\|C+\theta BB^{\intercal}\|)+\alpha}{\lambda_{min}(A)}$ with $\alpha>0$. If $\lambda_{1}=0$, we have $\lambda_{0}\leq-\lambda^{\prime}_{+}$ with $\lambda^{\prime}_{+}=\min(\lambda_{min}(C+BA^{-1}B^{\intercal}),\min(\|C\|,\|C
+\theta BB^{\intercal}\|)+\alpha)$.
- d.
Let $L=\max(\|A\|,\|B\|,\|C+\theta BB^{\intercal}\|)$. If $\tau\geq 1$, we have $\lambda_{0}^{2}+\lambda_{1}^{2}\leq\|M\|^{2}\leq 4\tau^{2}L^{2}$.

###### Proof.

Assume that $\lambda$ is an eigenvalue of $M$,
i.e. $\mbox{det}(\lambda I-M)=0$.
If $\lambda I+\tau A$ is invertible,
we can compute $\mbox{det}(\lambda I-M)$
using the Schur complement (54) to obtain

$$ $\displaystyle\mbox{det}(\lambda I-M)=\mbox{det}(\lambda I+\tau A)\mbox{det}(H( \lambda)).$ $$

Let the spectral decomposition of $A$ be $U\Lambda_{A}U^{\intercal}$, with orthogonal $U\in\mathbb{R}^{m\times m}$ and diagonal $\Lambda_{A}$. Then $H(\lambda)$ can be rewritten into (55).

Part a:
If $|\lambda_{1}|>\sqrt{\tau}\sqrt{1-\theta\lambda_{min}(A)}\|B\|$,
then $\lambda I+\tau A$ is invertible since $\lambda_{1}I$ is.
Therefore $\mbox{det}(H(\lambda))=0$. Consider

$$ $\displaystyle\frac{1}{\lambda_{1}}\mbox{Im}(H(\lambda))$ $\displaystyle=I+\frac{1}{\lambda_{1}}\tilde{B}\mbox{Im}(D)\tilde{B}^{\intercal}$ $\displaystyle=I-\frac{1}{\tau}\tilde{B}\mbox{diag}\left(\frac{1-\theta\lambda_ {k}(A)}{(\lambda_{0}/\tau+\lambda_{k}(A))^{2}+(\lambda_{1}/\tau)^{2}}\right)_{ k\leq m}\tilde{B}^{\intercal}$ $\displaystyle\geq I-\frac{1}{\tau}\frac{\tau^{2}}{\lambda_{1}^{2}}(1-\theta \lambda_{min}(A))\tilde{B}\tilde{B}^{\intercal}=I-\frac{\tau}{\lambda_{1}^{2}} (1-\theta\lambda_{min}(A))BB^{\intercal}$ $$

As $|\lambda_{1}|>\sqrt{\tau}\sqrt{1-\theta\lambda_{min}(A)}\|B\|$,
we have that $I-\frac{\tau}{\lambda_{1}^{2}}(1-\theta\lambda_{min}(A))BB^{\intercal}$ is p.d. and therefore
$\frac{1}{\lambda_{1}}\mbox{Im}(H(\lambda))$ is p.d., which leads to a contradiction.

Part b: Assume $\lambda_{0}>-\lambda_{+}$. As $\lambda_{1}\neq 0$, $\lambda I+\tau A$ is invertible
and $\mbox{det}(H(\lambda))=0$.
Following (57), we consider

$$ $\mbox{Re}(H(\lambda))+\frac{\tau\beta}{\lambda_{1}}\mbox{Im}(H(\lambda))=( \lambda_{0}+\tau\beta)I+C+\theta BB^{\intercal}+\tilde{B}F\tilde{B}^{\intercal}.$ $$

Take $\beta=\frac{1}{2\tau}(\lambda_{min}(A)\tau+\min(\|C\|,\|C+\theta BB^{\intercal
}\|))>0$,
then $\tilde{B}F\tilde{B}^{\intercal}$ is semi-p.d (positive semidefinite) because
each $f_{k}$ in (58) is larger than zero ($f_{k}\geq 0$).
Indeed, $\forall 1\leq k\leq m$,

$$ $\displaystyle\lambda_{0}/\tau+\lambda_{k}(A)-\beta$ $\displaystyle>-\lambda_{+}/\tau+\lambda_{k}(A)-\beta$ $\displaystyle=-\frac{\lambda_{min}(A)\tau-\min(\|C\|,\|C+\theta BB^{\intercal} \|)}{2\tau}-\beta+\lambda_{k}(A)$ $\displaystyle=\lambda_{k}(A)-\lambda_{min}(A)\geq 0.$ $$

Furthermore, either $(\lambda_{0}+\tau\beta)I+C$ or
$(\lambda_{0}+\tau\beta)I+C+\theta BB^{\intercal}$ is p.d. due
to the following:

- •
If $\|C\|<\|C+\theta BB^{\intercal}\|$,
we have $\lambda_{0}+\tau\beta>\|C\|$ because
$-\lambda_{+}+\tau\beta=\|C\|$.
- •
If $\|C\|\geq\|C+\theta BB^{\intercal}\|$,
we have $\lambda_{0}+\tau\beta>\|C+\theta BB^{\intercal}\|$ because
$-\lambda_{+}+\tau\beta=\|C+\theta BB^{\intercal}\|$.

As a consequence, $\mbox{Re}(H(\lambda))+\frac{\tau\beta}{\lambda_{1}}\mbox{Im}(H(\lambda))$ is p.d. which is a contradiction. Therefore $\lambda_{0}\leq-\lambda_{+}$.

Part c:
From Proposition J.1,
we know that $\lambda_{0}<0$.
If $\lambda_{1}=0$ and $0>\lambda_{0}>-\lambda^{\prime}_{+}$, $\lambda_{0}I+\tau A$ remains p.d. since

$$ $\lambda_{0}I+\tau A\geq(\tau\lambda_{min}(A)+\lambda_{0})I>(\min(\|C\|,\|C+ \theta BB^{\intercal}\|)+\alpha-\lambda^{\prime}_{+})I\quad\mbox{semi-p.d.}$ $$

Therefore $\lambda I+\tau A$ is invertible
and $\mbox{det}(H(\lambda))=0$.
Following (55), we have

$$ $\displaystyle H(\lambda)$ $\displaystyle=\lambda_{0}I+C+\theta BB^{\intercal}+\tilde{B}\mbox{diag}\left( \frac{1-\theta\lambda_{k}(A)}{\lambda_{0}/\tau+\lambda_{k}(A)}\right)_{k\leq m }\tilde{B}^{\intercal}$ $\displaystyle=\lambda_{0}I+C+\theta BB^{\intercal}+BA^{-1}B^{\intercal}-\tilde {B}\mbox{diag}\left(\frac{\lambda_{0}/\tau+\theta\lambda_{k}^{2}(A)}{(\lambda_ {0}/\tau+\lambda_{k}(A))\lambda_{k}(A)}\right)_{k\leq m}\tilde{B}^{\intercal}$ (60) $\displaystyle=\lambda_{0}I+C+BA^{-1}B^{\intercal}-\tilde{B}\mbox{diag}\left( \frac{(1-\theta\lambda_{k}(A))\lambda_{0}/\tau}{(\lambda_{0}/\tau+\lambda_{k}( A))\lambda_{k}(A)}\right)_{k\leq m}\tilde{B}^{\intercal}$ (61) $\displaystyle>C+BA^{-1}B^{\intercal}-\lambda^{\prime}_{+}I$ $\displaystyle\geq C+BA^{-1}B^{\intercal}-\lambda_{min}(C+BA^{-1}B^{\intercal}) I\quad\mbox{semi-p.d.}$ $$

As a consequence, $H(\lambda)$ is p.d. which is a contradiction. Therefore $\lambda_{0}\leq-\lambda^{\prime}_{+}$.

Part d: Note that $\|I-\theta A\|\leq 1$.
For any $(x,y)\in\mathbb{R}^{n}\times\mathbb{R}^{m}$, we compute

$$ $\displaystyle\left\|M\left(\begin{array}[]{l}x\\ y\end{array}\right)\right\|$ $\displaystyle=\left\|\begin{array}[]{l}-Cx-By-\theta BB^{\intercal}x+\theta BA ^{\intercal}y\\ \tau B^{\intercal}x-\tau Ay\end{array}\right\|$ $\displaystyle=\sqrt{\|-Cx-By-\theta BB^{\intercal}x+\theta BA^{\intercal}y\|^{ 2}+\|\tau B^{\intercal}x-\tau Ay\|^{2}}$ $\displaystyle\leq\sqrt{2(\|(C+\theta BB^{\intercal})x\|^{2}+\|(B-\theta BA^{ \intercal})y\|^{2})+2\tau^{2}(\|B^{\intercal}x\|^{2}+\|Ay\|^{2})}$ $\displaystyle\leq\sqrt{2\|C+\theta BB^{\intercal}\|^{2}\|x\|^{2}+2\|B\|^{2}\|I -\theta A\|^{2}\|y\|^{2}+2\tau^{2}\|B^{\intercal}\|^{2}\|x\|^{2}+2\tau^{2}\|A \|^{2}\|y\|^{2}}$ $\displaystyle\leq\sqrt{2(1+\tau^{2})L^{2}(\|x\|^{2}+\|y\|^{2})}\leq 2\tau L \left\|\left(\begin{array}[]{l}x\\ y\end{array}\right)\right\|.$ $$

Therefore $\|M\|\leq 2\tau L$ and $\lambda_{0}^{2}+\lambda_{1}^{2}\leq\|M\|^{2}\leq 4\tau^{2}L^{2}$.
∎

From Lemma J.1,
we are ready to obtain an upper bound
of the local convergence rate of Asymptotic $\tau$-SGA.
Recall that $L_{s}=\max(\|\mathbf{A}\|,\|\mathbf{B}\|,\|\mathbf{C}+\theta\mathbf{B}\mathbf{
B}^{\intercal}\|)$
and $\mu_{s}=\min(L_{s},\lambda_{min}(\mathbf{C}+\mathbf{B}\mathbf{A}^{-1}\mathbf{B
}^{\intercal}))$.

###### Proposition J.2 .

Assume $(x^{\ast},y^{\ast})$ is a DSE of $f\in C^{2}$.
If $\tau\geq\frac{2L_{s}}{\lambda_{min}(\mathbf{A})}$
and $\gamma=\frac{1}{4\tau L_{s}}$, we have
$\rho(I+\gamma M_{s})\leq 1-\frac{\mu_{s}}{16\tau L_{s}}$.

###### Proof.

Let $L=\max(\|\bar{\mathbf{A}}\|,\|\bar{\mathbf{B}}\|,\|\bar{\mathbf{C}}+\theta\bar
{\mathbf{B}}\bar{\mathbf{B}}^{\intercal}\|)$
and $\mu_{x}=\min(L,\lambda_{min}(\bar{\mathbf{C}}+\bar{\mathbf{B}}\bar{\mathbf{A}}
^{-1}\bar{\mathbf{B}}^{\intercal}))$.
We verify that $\lambda_{min}(\mathbf{A})=\lambda_{\min}(\bar{\mathbf{A}})$, $L=L_{s}$ and $\mu_{x}=\mu_{s}$.
Let $\lambda=\lambda_{0}+i\lambda_{1}$ be an eigenvalue of $M_{s}$.
It follows that we only need to show that $|1+\gamma\lambda|\leq 1-\frac{\mu_{x}}{16\tau L}$.

To apply Lemma J.1, we denote $A=\bar{\mathbf{A}}$, $B=\bar{\mathbf{B}}$ and $C=\bar{\mathbf{C}}$. We rewrite that
$\tau\lambda_{min}(A)=\min(\|C\|,\|C+\theta BB^{\intercal}\|)+\alpha$
with $\alpha\geq 2L-\min(\|C\|,\|C+\theta BB^{\intercal}\|)>0$.
There are two cases to verify:

- •
Case $\lambda_{1}=0$:
we have $\lambda_{0}\leq-\lambda_{+}^{\prime}$. Since $\tau\lambda_{min}(A)>L$, we have
$\lambda_{+}^{\prime}=\min(\tau\lambda_{min}(A),\lambda_{min}(C+BA^{-1}B^{
\intercal}))\geq\min(L,\lambda_{min}(C+BA^{-1}B^{\intercal}))=\mu_{x}$.
Thus $\lambda_{0}\leq-\mu_{x}$ and $1+\gamma\lambda_{0}\leq 1-\frac{1}{4\tau L}{\mu_{x}}$.
On the other hand, $\tau\geq 2$ and $\lambda_{0}\geq-2\tau L$, thus $1+\gamma\lambda_{0}\geq 1-\frac{1}{4\tau L}{2\tau L}=1/2$. Therefore $|1+\gamma\lambda_{0}|\leq 1-\frac{1}{4\tau L}{\mu_{x}}\leq 1-\frac{\mu_{x}}{16
\tau L}$.
- •
Case $\lambda_{1}\neq 0$: we know that $-2\tau L\leq\lambda_{0}\leq-\lambda_{+}=-\frac{\alpha}{2}$ and $|\lambda_{1}|\leq\sqrt{\tau}L$. Note that $\alpha=\tau\lambda_{min}(A)-\min(\|C\|,\|C+\theta BB^{\intercal}\|)\geq\tau
\lambda_{min}(A)-L$ and by assumption $L\leq\frac{\tau\lambda_{min}(A)}{2}$. Thus $\alpha\geq\frac{\tau\lambda_{min}(A)}{2}$. It follows that
$\displaystyle|1+\gamma\lambda|^{2}$
$\displaystyle=(1+\gamma\lambda_{0})^{2}+\gamma^{2}\lambda_{1}^{2}$
$\displaystyle\leq\left(1-\gamma\frac{\alpha}{2}\right)^{2}+\gamma^{2}\tau L^{2}$
$\displaystyle=\left(1-\frac{\alpha}{8\tau L}\right)^{2}+\frac{1}{16\tau}$
$\displaystyle\leq\left(1-\frac{\lambda_{min}(A)}{16L}\right)^{2}+\frac{\lambda
_{min}(A)}{32L}\leq 1-\frac{\lambda_{min}(A)}{16L}.$
Therefore $|1+\gamma\lambda|\leq\sqrt{1-\frac{\lambda_{min}(A)}{16L}}\leq 1-\frac{\lambda
_{min}(A)}{32L}\leq 1-\frac{L}{16\tau L}\leq 1-\frac{\mu_{x}}{16\tau L}$.

∎

###### Lemma J.1 .

###### Proof.

###### Proposition J.2 .

###### Proof.

## Appendix K Computational efficiency of τ 𝜏 \tau italic_τ -SGA

We first discuss the computation of $\tilde{\mathbf{B}}[\eta]$ and $\tilde{\mathbf{B}}^{\intercal}[\delta]$
when $\mathcal{M}_{1}$ (resp. $\mathcal{M}_{2}$) is an embedded
sub-manifold of $\mathbb{R}^{d_{1}^{\prime}}$ (resp. $\mathbb{R}^{d_{2}^{\prime}}$).
We then propose a linear-time computational procedure using auto-differentiation when $\mathcal{M}_{1}$ is Euclidean,
which is applicable to Wasserstein GAN.
When $\mathcal{M}_{2}$ is also Euclidean, this
procedure is equivalent to the one proposed in Balduzzi et al. ().

To compute $\tilde{\mathbf{B}}[\eta]$ at $(x,y)$,
we first fix $\eta=\eta(x,y)\in T_{y}\mathcal{M}_{2}$.
From the property of embedded
sub-manifolds, we use $\leftrightarrow$ to identify a tangent vector
$\delta\in T_{x}\mathcal{M}_{1}$
with a vector $\delta^{\prime}=(\delta^{\prime}_{i})_{i\leq d_{1}^{\prime}}\in\mathbb{R}^{d_{
1}^{\prime}}$
(resp. $\eta\in T_{y}\mathcal{M}_{2}$ with $\eta^{\prime}=(\eta^{\prime}_{i})_{i\leq d_{2}^{\prime}}\in\mathbb{R}^{d_{2}^{
\prime}}$).

This allows one to compute the cross-gradients in the embedded space, by

$$ $\tilde{\mathbf{B}}[\eta](x,y)=D_{y}\mbox{grad}_{x}f(x,y)[\eta]\quad \leftrightarrow\quad D_{y}\delta^{\prime}(x,y)[\eta^{\prime}]=(\langle(\mbox{ grad}_{y}\delta^{\prime}_{i}(x,y))^{\prime},\eta^{\prime}\rangle_{y})_{i\leq d _{1}^{\prime}}\in\mathbb{R}^{d_{1}^{\prime}}.$ $$

Note that $\langle(\mbox{grad}_{y}\delta^{\prime}_{i}(x,y))^{\prime},\eta^{\prime}\rangle
_{y}$
is computed on $\mathbb{R}^{d_{2}^{\prime}}$
with a metric induced from $T_{y}\mathcal{M}_{2}$.
Assume that its computational time
is $O(d_{2}^{\prime})$ for each $i$.
Then it takes $O(d_{2}^{\prime}d_{1}^{\prime})$ to compute $\tilde{\mathbf{B}}[\eta](x,y)$ in the embedded space.
We can obtain a similar cost for $\tilde{\mathbf{B}}^{\intercal}[\delta](x,y)$.

When $\mathcal{M}_{1}=\mathbb{R}^{d_{1}}$, the computational complexity of $\tilde{\mathbf{B}}[\eta](x,y)$
can be significantly reduced: $\forall i\leq d_{1}$,

$$ $D_{y}\partial_{x_{i}}f(x,y)[\eta]=\partial_{x_{i}}\langle\mbox{grad}_{y}f(x,y) ,\eta\rangle_{y}\quad\leftrightarrow\quad\partial_{x_{i}}\langle\eta^{\prime}( x,y),\eta^{\prime}\rangle_{y}.$ (62) $$

Importantly, one does not need to recompute $\langle\eta^{\prime}(x,y),\eta^{\prime}\rangle_{y}$ for each $i$.
Therefore, the whole cost of $\tilde{\mathbf{B}}[\eta](x,y)$ is $O(d_{1}+d_{2}^{\prime})$.
Note that in (62), the term $\eta^{\prime}$ is “detached”.

For $\tilde{\mathbf{B}}^{\intercal}[\delta](x,y)$, we “detach” $\delta=\delta(x,y)\in T_{x}\mathcal{M}_{1}$
and compute

$$ $\tilde{\mathbf{B}}^{\intercal}[\xi_{1}](x,y)=\sum_{i\leq d_{1}}\partial_{x_{i} }\eta(x,y)\delta_{i}\quad\leftrightarrow\quad\partial_{\eta^{\prime\prime}} \left(\sum_{i\leq d_{1}}\partial_{x_{i}}\langle\eta^{\prime}(x,y),\eta^{\prime \prime}\rangle\delta_{i}\right)_{|_{\eta^{\prime\prime}=0}}.$ $$

This implies that we can first compute $\partial_{x_{i}}\langle\eta^{\prime}(x,y),\eta^{\prime\prime}\rangle$ for each $i$ as
in (62). We then compute its sum with $\delta_{i}$ which takes $O(d_{1})$.
Finally an extra auto-differentiation is taken with respect to $\eta^{\prime\prime}$ which costs $O(d_{2}^{\prime})$.
The whole cost of $\tilde{\mathbf{B}}^{\intercal}[\xi_{1}](x,y)$ is therefore $O(d_{1}+d_{2}^{\prime})$.
Note that in this case, we use the Euclidean metric on $\mathbb{R}^{d_{2}^{\prime}}$ to evaluate $\langle\eta^{\prime}(x,y),\eta^{\prime\prime}\rangle$ rather than the induced Riemannian metric in (62).

### K.1 Extension to stochastic τ 𝜏 \tau italic_τ -SGA

We construct stochastic $\tau$-SGA
through an unbiased estimation of the terms in the update rule
of deterministic $\tau$-SGA.
To achieve this,
we compute
$\partial_{x_{i}}\langle\eta^{\prime}(x,y),\eta^{\prime}\rangle_{y}$
in (62)
using two mini-batches independently sampled from a training set,
one to estimate $\eta^{\prime}(x,y)$, the other
to estimate $\eta^{\prime}$.
Similarly,
we use these mini-batches
to estimate
$\eta^{\prime}(x,y)$
and $\delta$ in
$\tilde{\mathbf{B}}^{\intercal}[\xi_{1}](x,y)$.

### K.2 Numerical results on computational time

In Figure 4,
we compare their computation time
in terms of the value of $f$,
or the FID (test) score
in Example 2 of Section 2.2
(with deterministic gradients)
and the image datasets in Section 4
(with stochastic gradients).

In Example 2,
we find that $\tau$-SGA has a significant speedup compared to $\tau$-GDA.
On the two image datasets, the speedup is less significant. On the MNIST (digit 0) dataset, we find that $\tau$-SGA is slightly faster and it can reach a lower FID (test) score. On Fashion-MNIST, the speed of $\tau$-GDA and $\tau$-SGA is similar.
In these cases, we find that the computational time per iteration (per $t$)
of $\tau$-SGA is roughly 3-4 times that of $\tau$-GDA.
A further improvement of the computational efficiency of
$\tau$-SGA can make its speedup more significant.

Figure: (a)
Refer to caption: extracted/5893895/plots/neurips_compute_dse_dert_sgatau_speedup_example2.png

## Appendix L Details of numerical experiments

### L.1 Choice of smooth non-linearity

We aim to build GAN models
whose value function $f$ is
twice continuously differentiable,
based on smooth non-linearities studied in Biswas et al. ().
For the discriminator of Gaussian distribution,
$\rho$ is a smooth approximation
of the absolute value non-linearity of the form

$$ $\sigma(a)=\sqrt{a^{2}+\epsilon^{2}},\quad\epsilon=10^{-6}.$ $$

This function is twice continuously differentiable
since $\sigma^{\prime\prime}(a)=\frac{\epsilon^{2}}{2(a^{2}+\epsilon^{2})^{3/2}}$.
Similarly, for the discriminator of MNIST digits,
$\sigma$ is a smooth ReLu non-linearity,

$$ $\sigma(a)=(a+\sqrt{a^{2}+\epsilon^{2}})/2.$ $$

### L.2 DCGAN generator for MNIST

We use the default DCGAN generator to model MNIST (digit 0),
implemented in WGAN-GP(^6^66[https://github.com/caogang/wgan-gp](https://github.com/caogang/wgan-gp)).
We choose this dataset because there is no batch normalization module
used in the generator network (which is typically used on other datasets such as CIFAR-10).
From the 60k training samples of MNIST,
we take $4932$ training samples, $991$ validation samples
of $\phi_{data}$ to build the dataset MNIST (digit 0).
We also make a slight modification of the generator
so that the function $x\mapsto G_{x}(Z)$
is twice continuously differentiable at any $x\in\mathcal{M}_{1}$
for a fixed $Z$.
For this, each ReLu non-linearity in $Z\mapsto G_{x}(Z)$ is replaced
by the smooth ReLu non-linearity
in Appendix L.1.

### L.3 Scattering CNN discriminator for MNIST

We construct a smooth Lipschitz-continuous discriminator
for MNIST digits,

$$ $D_{y}(\phi)=\langle v_{y},\sigma(w_{y}\star P(\phi)+b_{y})\rangle,$ $$

where $\rho$ is the smooth ReLu defined in Appendix L.1.
For a fixed $\phi$, this makes the function $y\mapsto D_{y}(\phi)$
twice continuously differentiable at any $y\in\mathcal{M}_{2}$.
However, the function $\phi\mapsto D_{y}(\phi)$
is not everywhere twice continuously differentiable
due to the modulus non-linearity in the scattering transform.
We therefore replace this modulus non-linearity $z\mapsto|z|=|z_{re}+iz_{im}|$
by $z\mapsto\sqrt{z_{re}^{2}+z_{im}^{2}+\epsilon^{2}}$
for each complex number input $z=z_{re}+iz_{im}$.
This makes the function $\phi\mapsto D_{y}(\phi)$ twice continuously differentiable.
As a consequence, the function $(x,y)\mapsto f(x,y)$
is also twice continuously differentiable,
which is induced from the smoothness
of $(x,y)\mapsto D_{y}(G_{x}(Z))$
and $y\mapsto D_{y}(\phi_{data})$.

##### Scattering transform

The input $\phi$ with dimension $d=784$ is represented as an image
of size $28\times 28$.
It is pre-processed by the wavelet scattering transform $P(\phi)$
to extract stable edge-like information
using Morlet wavelet of different orientations and scales.
We use the second-order scattering transform
with four wavelet orientations (between $[0,\pi)$)
and two wavelet scales.
It first computes
the convolution of $\phi$ with
each wavelet filter,
then a smooth modulus non-linearity is applied
to each feature map.
This computation is repeated one more time
on each obtained feature maps
and then a low-pass filter is applied to
each of the channels.
The obtained scattering features $P(\phi)$
is an image of size $9\times 9$ with $25$ channels ($I=25$,$n=9$).

##### Orthogonal CNN layer

The orthogonal CNN layer is parameterized by
the kernel $w_{y}$ and bias $b_{y}$.
The kernel $w_{y}$ has $5\times 5$ spatial size ($k=5$).
With a suitable padding and stride (two by two), we obtain
an output image of size $5\times 5$ ($N=5$) with $J=256$ channels.
Therefore the embedding space dimension
of $v_{y}$ is $JN^{2}=6400$.

### L.4 Stiefel manifold geometry

For the discriminators of the Gaussian distribution and MNIST (digit 0),
part of the parameters in $y$ belong to Stiefel manifolds.
To choose a Riemannian metric on a Stiefel manifold (which is non-Euclidean),
we use the one in Manton (, equation 20).
We also use the SVD projection in Manton (, Proposition 12)
as the retraction $\mathcal{R}_{2}$ on each Stiefel manifold.

### L.5 Initialization for local convergence

##### Simultaneous τ 𝜏 \tau italic_τ -GDA initialization for Gaussian distribution

Starting from a random initialization of $(x,y)$,
we apply the stochastic $\tau$-GDA method to build a pre-trained model.
It is pre-trained with batch size $1000$, learning rate $\gamma=0.0002$ and $\tau=100$
for $T=50000$ iterations.

##### Alternating τ 𝜏 \tau italic_τ -GDA initialization for MNIST (digit 0)

The stochastic alternating $\tau$-GDA (Zhang et al., )
is often used in the training of WGAN-GP (Gulrajani et al., ),
and it can be extended to Riemannian manifold using a suitable retraction.
Starting from a random initialization of $(x,y)$,
we apply the alternating $\tau$-GDA method to build a pre-trained model.
Each iteration amounts to
perform $\tau=5$ gradient updates of $y$
(with learning rate $0.1$) and one gradient
update of $x$ (with learning rate $0.1$).
It is pre-trained with batch size $128$
for a total number of $10000$ iterations.
We obtain a pre-trained model
with FID (train) $=12$ and FID (val) $=18$.

### L.6 Extra numerical experiments

To make the results in Section 4 more convincing,
we add one experiment on a subset of the Fashion-MNIST dataset (restricted to the T-shirt class).
We use the same GAN architecture from the MNIST (digit 0) dataset.

From the results in Table 3, we see that the behavior of the $\tau$-GDA and $\tau$-SGA algorithms (in terms of the f and angle measures) are similar to the MNIST (digit 0) case. As in Table 2, we observe that one can improve the local convergence of $\tau$-GDA in the range of small $\tau$ by using $\tau$-SGA. We also find that an improved convergence (in terms of the angle)
lead to improved GAN models in terms of the FID scores.

**Table 3: Last iteration measures of stochastic $\tau$-GDA and $\tau$-SGA on the Wasserstein GAN of Fashion-MNIST (T-shirt). We report the $f$, angle, and FID scores computed at $(x(T),y(T))$. $\tau$-GDA is trained for $T=2\times 10^{5}$ iterations with $\gamma=0.1/\tau$. $\tau$-SGA is trained longer with $\gamma=0.02$, $\tau=5$, $\theta=0.075$. Only significant digits are reported with respect to the standard deviation (shown in $(\pm)$) of the statistical estimation of these quantities.**
|  | $\tau$-GDA |  |  |  |
| --- | --- | --- | --- | --- |
| $\tau$ | $f$ | angle | FID (train) | FID (val) |
| $5$ | -0.52 | -0.31 | 80 ($\pm 0.4$) | 92 ($\pm 1$) |
| $10$ | 0.03 ($\pm 0.03$) | 0.05 ($\pm 0.04$) | 73 | 86 |
| $30$ | 0.01 | 0.15 ($\pm 0.03$) | 24.7 ($\pm 0.2$) | 38.7 ($\pm 0.5$) |

## Peer Reviews

confidence: 3
contribution: 3
initial_rating: 6
presentation: 2
questions: - lines 45-47 could be improved to better explain this paper's focus on differential equilibriums (DNE and DSE). This could have been local Nash and local minimax points. Then why differential equilibrium? 
- Table 1 is confusing: What is your focus here? The title should be more specific to make a point. Isn't Nash equilibrium found by Zhang et al. (2023) DNE? This is not clear in the table. Also, your result is local (both in optimality and convergence) while others might be global. I am not sure whether local/global should be discussed in the table, but I am sure that this table needs to be reorganized and rewritten.
- I suggest rewriting Definition 3.1 to make it more explicit that it defines local linear convergence.
- line 276: Were you trying to say that $M$ with eigenvalues with only negative real parts?
- lines 524-525: I don't follow what do you mean here.
soundness: 3
strengths: Existing results on local convergence to DSE and DNE in nonconvex-nonconcave min-max problems in Euclidean space are crucial. Thus, while its extension to Riemannian manifolds may seem straightforward, it is a significant and non-trivial advancement. Although I was not able to check all the details of the proof, the parts I checked looked correct.
summary: This paper extends the analysis of the local convergence of two-time-scale $\tau$-GDA and $\tau$-SGA methods to the differential Stakelberg equilibrium (DSE) of Euclidean min-max problem to the Riemannian min-max problem. This paper first generalizes the notion of DSE in the Euclidean space to the Riemannian manifold (in Section 2.1), using a local coordinate chart. Then, the authors provide a sufficient condition (in Theorem 3.1) for the fixed point method to have a local linear convergence. Based on Theorem 3.1, this paper shows a specific sufficient condition for the $\tau$-GDA to locally linearly converge to DSE in Theorem 3.2. In addition, the authors generalize $\tau$-SGA for Riemannian min-max problem, which is originally developed to mitigate undesirable rotational dynamics encountered by $\tau$-GDA in the Euclidean space. In the Euclidean space, extragradient (or optimistic gradient) method is widely used to relieve such rotational behavior, but its Riemannian version, named RCEG, is computationally expensive. This paper claims that their $\tau$-SGA equipped with auto-differentiation may have computational benefit over RCEG. Experiments on a toy example and a more realistic Orthogonal WGAN problem empirically confirm this paper's theoretical findings.
weaknesses: Missing preliminaries for non-experts on Riemannian manifold: For example, adding the definition of Riemannian gradient/Hessian/cross-gradient, the role of local coordinate chart and the relation between $f$ and $\bar{f}$ would help readers to better understand the context.

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 4
questions: From the experiments, it seems like τ-SGA is used as a sort of “last mile” algorithm to guarantee fast convergence after good initialization/pretraining (either via ansatz in the analytical examples of Figure 1 or via τ-GDA for the GANs). Is this a reasonable interpretation of the authors’ proposed use of τ-SGA, and is there any reason that τ-GDA—in particular, alternating τ-GDA—is used to pretrain the GANs?

In Theoerems 3.2 and 3.4, after the main statements there are additional (“furthermore…”) statements for when γ is fixed exactly depending on L_g/L_s. Is this reasonable given that L_g depends on the DSE, which is not known in advance? (Maybe this is not so important, but as it is written I am not sure what the significance of these “furthermore” statements is.)

τ-SGA has an additional hyperparameter θ (or μ) to tune compared to τ-GDA, and in the GAN example of Table 2 a “suitable choice” of θ seems to have been made. Does the value of θ matter much, and how is it chosen?

I am not so familiar with the way “retraction” is used in Section 3.1, though it does not affect much as the algorithms are easily understood. Perhaps some quick definition or source could be given?
soundness: 4
strengths: Overall, this paper seems like a significant algorithmic contribution with sound theoretical grounding. The authors extend known, tight local convergence results for τ-GDA from Euclidean space to Riemannian manifolds and introduce τ-SGA to the context of Riemannian manifolds. Both in theory and practice τ-SGA is shown convincingly to be able to achieve (local convergence) in settings where τ-GDA cannot, highlighting its usefulness in this setting of Riemannian min-max problems with minimal assumptions (just f twice continuously differentiable).

The paper is generally well written. The sections on differentiable equilibria, algorithms, and local convergence read especially well, with much background and explanation given to readers. Motivations and the context of prior work and clearly written.
summary: This paper analyzes the local convergence of the deterministic τ-GDA and τ-SGA algorithms to differentiable equilibria in min-max games on Riemannian manifolds. Using the Ostrowski theorem, sufficient conditions on the algorithm hyperparameters (in particular, on the learning ratio τ) are given for linear-rate local convergence to Differentiable Stackelberg Equilibria (DSE) for both τ-GDA and Asymptotic τ-SGA (an approximation to τ-SGA that is amenable to easier analysis). Local convergence to Differentiable Nash Equilibria (DNE) is also given for τ-GDA.

This paper marks the first time that τ-SGA has been written for Riemannian manifolds, and the authors put it forth as an algorithm that can solve the problem of τ-GDA’s slow convergence rate in certain settings. Both theory and example are given to show that Asymptotic τ-SGA can be locally convergent at a smaller τ compared to τ-GDA, and with a faster rate; moreover, in Euclidean space the rate guarantee outperforms that of the extra-gradient method, which has been used to overcome the slow convergence rate of τ-GDA but is costly in Riemannian space. Experiments also show that the gap between Asymptotic τ-SGA and the true τ-SGA is small, so the theoretical results should also hold for the more natural τ-SGA.

Finally, the paper applies stochastic τ-GDA and τ-SGA to train orthogonal Wasserstein GANs with a discriminator parametrized in the Stiefel manifold. It is shown that even with good initialization a small τ may cause τ-GDA to oscillate with high amplitude, but that τ-SGA can converge even at those small values of τ.
weaknesses: While it is clearly stated in the introduction that global convergence is difficult to achieve in this setting and local convergence (to differentiable equilibria) may be all we can get theoretically, it is not clear what these limitations mean, or are expected to mean, in practice. The experiments are set up to show local convergence after good initialization, but it is not clear what the good initialization entails or if τ-SGA can be used when you are not near a DSE. Also maybe some comment could be made about the existence and significance of DSE in the settings provided (e.g., the Wasserstein GANs), and about other methods for doing min-max Riemannian optimization in these setting (e.g., over Steifel manifolds).

I also find the writing of the GAN experiments and their setups/parameters/initialization a bit terse and less clearly written than what comes before, with a few things in particular lacking explanation. I ask some questions about this below.

---

confidence: 2
contribution: 3
initial_rating: 6
presentation: 3
questions: Question for Theorem 3.1. (Ostrowski Theorem on manifold)

- The author indicates that the theorem is analogous to the theorem in the Euclidean case. Is the Euclidean version of the theorem also based on the spectral radius of a linear transformation similar to (9)? 

- What is the size of the local convergence radius around ( 𝑥*, 𝑦*)? 

- How can one ensure that the initialization lies within this local neighborhood? Does the initialization used in experiments satisfy this requirement? 

- What is the reasoning behind assuming a small spectral radius for 𝐷𝑇* ?
soundness: 3
strengths: This article provides sufficient conditions for the local convergence of τ-GDA and τ-SGA to differential equilibrium on a Riemannian manifold, where the author proves a linear convergence rate of τ-GDA to DSE and DNE, depending on the spectral radius of the Jacobian at equilibrium. THe author introduces τ-SGA on Riemannian manifolds to further improve this rate.
summary: Using the concepts of DSE and DNE on Riemannian manifolds, the author derives conditions on the range of τ and learning rate of x to ensure linear convergence of τ-GDA to DSE and DNE.

It introduces a novel algorithm τ-SGA to enhance the convergence of τ-GDA, allowing a broader range of τ for convergence when τ is large, and achieving faster local convergence than τ-GDA.

Applying these insights, the author improves the training of orthogonal Wasserstein GANs. Numerical results show that enhanced convergence of stochastic τ-GDA and τ-SGA can improve the learned generator on simple benchmarks.
weaknesses: - Some key terms could benefit from further explanation or intuitive insights. For example, what are the underlying intuitions behind the linear transformations in (9) and (10)?

- It appears that several theorems are extensions of those from the Euclidean case. What are the main differences or key insights that allowed these results to be successfully extended to the Riemannian manifold?

- Given the author's claim of linear convergence rates for the proposed methods, are there any experimental results to support this?

---

confidence: 1
contribution: 1
initial_rating: 1
presentation: 2
questions: I think that this manuscript contributes nothing in science expect some hard understanding of mathematical symbols.
soundness: 2
strengths: No
summary: This paper studies min-max algorithms to solve zero-sum differential games on Riemannian
manifold
weaknesses: This paper does not have any interesting parts. Firstly, the Riemannian optimization is not practical. Furthermore, could you show some practical utilization of your algorithms?

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 3
questions: 1. In line 277, what is the order of $\gamma^{\cdot}(M)$? Could the authors provide an intuitive explanation of this function?
2. The computation of $\tau$-SGA involves a matrix-vector product. Are the benefits sufficient to justify the additional computational cost?
soundness: 3
strengths: 1. **Significance**: This paper addresses the more general nonconvex-nonconcave setting, extending previous results to Riemannian manifolds.
2. **Comprehensiveness**: It provides a range of results by examining two types of equilibrium points, two algorithms, and extensive numerical results.
summary: This paper investigates the local convergence of simultaneous minimax algorithms to differential Stackelberg equilibrium (DSE) and differential Nash equilibrium (DNE) on Riemannian manifolds. The authors study two representative deterministic algorithms, $\tau$-GDA and $\tau$-SGA. For $\tau$-GDA, they derive sufficient conditions for convergence to DSE and DNE. For $\tau$-SGA, they analyze an asymptotic variant and demonstrate a faster convergence rate to DSE under specific conditions. Furthermore, the authors conduct numerical experiments to explore the behavior of stochastic versions of the two algorithms.
weaknesses: 1. **Limited Comparative Analysis and Interpretation**: The relationship and comparative insights between the theoretical results could be strengthened. Here are some examples.
* (i) The paper examines two equilibrium points, DSE and DNE, where DNE is a subclass of DSE with stricter definitions. Convergence to DNE would be expected to require fewer restrictions or achieve a faster rate. Theorems 3.2 and 3.3 support the former expectation, as the range of $\tau$ in Theorem 3.3 is broader than in Theorem 3.2, yet no convergence rate is provided for Theorem 3.3.
* (ii) The authors seem to aim to show the superiority of $\tau$-SGA over  $\tau$-GDA. However, by comparing Theorems 3.2 and 3.4, the range of $\tau$ ensuring convergence for asymptotic $\tau$-SGA exceeds that of $\tau$-GDA only if $\| C + \theta BB^\top \| < \| C \|$. Further discussion on practical ways to ensure this condition or examples where it holds would be helpful. Additionally, there is no theoretical counterpart to Theorem 3.3 for $\tau$-SGA. While $\tau$-SGA is proposed to avoid rotational dynamics, this advantage is demonstrated only in numerical results (e.g., Figure 1), not theoretically, leaving theoretical support insufficient.

2. **Unclear Analysis for (Asymptotic) $\tau$-SGA**: 
* (i) In line 327, the asymptotic analysis requires $\theta$ to be of constant order. Why, then, are parameters $\mu$ and $\tau+1$ introduced in Eqs. 12 and 13, and what is their significance?
* (ii) If the authors aim to show the superiority of $\tau$-SGA, Theorem 3.4 alone is insufficient, and Figure 1's simple example is not enough to establish the similarity between $\tau$-SGA and its asymptotic variant. The claim "Theorem 3.4 is valid for $\tau$-SGA" in line 363 lacks rigor. If deriving a theoretical result for $\tau$-SGA is challenging, additional numerical experiments for the asymptotic $\tau$-SGA$ should at least be provided.

Meta-review: This paper extends the analysis of local convergence for deterministic $\tau$-GDA and $\tau$-SGA algorithms to min-max problems on Riemannian manifolds, generalizing concepts like Differentiable Stackelberg Equilibrium (DSE) and Differentiable Nash Equilibrium (DNE) from Euclidean spaces. Using theoretical tools like Ostrowski's theorem, the authors provide sufficient conditions for local linear convergence of $\tau$-GDA and a faster convergence rate for an asymptotic variant of $\tau$-SGA, which addresses the slow convergence and undesirable rotational dynamics of $\tau$-GDA. Unlike the computationally expensive Riemannian extragradient methods, $\tau$-SGA offers practical advantages with minimal performance gaps between its theoretical and natural forms. Empirical validation includes training orthogonal Wasserstein GANs on the Stiefel manifold, demonstrating $\tau$-SGA's robustness and improved convergence even under challenging parameter settings.

This paper seems like a significant algorithmic contribution with sound theoretical grounding. The authors extend known, tight local convergence results for τ-GDA from Euclidean space to Riemannian manifolds and introduce τ-SGA to the context of Riemannian manifolds. Existing results on local convergence to DSE and DNE in nonconvex-nonconcave min-max problems in Euclidean space are crucial. Thus, while its extension to Riemannian manifolds may seem straightforward, it is a significant and non-trivial advancement. Both in theory and practice τ-SGA is shown convincingly to be able to achieve (local convergence) in settings where τ-GDA cannot, highlighting its usefulness in this setting of Riemannian min-max problems with minimal assumptions (just f twice continuously differentiable).

The paper is generally well written. The sections on differentiable equilibria, algorithms, and local convergence read especially well, with much background and explanation given to readers. Motivations and the context of prior work and clearly written.

Decision: Accept (Poster)