Nonmonotone feasible arc search algorithm for minimization on Stiefel manifold

Abstract

We devise a new numerical method for solving the minimization problem over the Stiefel manifold, that is, the set of matrices of order $n \times p$ (here $p \le n$) with orthonormal columns. Our approach consists in a nonmonotone feasible arc search along a sufficient descent direction to assure convergence to stationary points, regardless the initial point considered. The feasibility of the iterates is maintained through a variation of the Cayley transform and thus our scheme can be seen as a retraction-based algorithm for minimization with orthogonality constraints. We emphasize that our scheme solves a $p\times p$ linear system at each iteration and has computational complexity of $O(np^2) + O(p^3)$, which is interesting when $p \ll n$. We present a general algorithmic framework for minimization on Stiefel manifold, give its global convergence properties and report numerical experiments on interesting applications.

Appendix A: Floating point arithmetic considerations

We remark that in Theorem 1 the calculations of ${{\mathcal {R}}}_T$, U and V are assumed to be computed with high precision, as well as the solution of the linear system (A2), in such a way that $X^\textrm{T}X = I$ throughout the proof. Indeed, this assumption is even observed in our numerical experiments. Nevertheless, occasionally, such precision can be lost and such a consideration is no longer true. In this direction, we also provide a result analogous to Theorem 1 but considering the feasibility residual ${{{\mathcal {R}}}}_V = I - X^\textrm{T}X$ into calculations. Therefore, in this context of loss of numerical precision, we recommend to replace ${{\mathcal {R}}}_1$ and ${{\mathcal {R}}}_T$ in (20) by ${{\mathcal {R}}}_1^\textrm{Full}$ and ${{\mathcal {R}}}_T^\textrm{Full}$, respectively, according to the next theorem.

Theorem 3

Let $X \in \Gamma$, $Z \in {\mathbb {R}}^{n\times p}$ and define $W(t)\in {\mathbb {R}}^{2p \times p}$, ${{\mathcal {R}}}_1$ and ${{\mathcal {R}}}_T$ as in Theorem 1. For all $t\in {\mathbb {R}}$, if we consider ${{\mathcal {R}}}_V = I - X^\textrm{T}X$,

$$\begin{aligned} {{\mathcal {R}}}_1^\textrm{Full} = \left( {{\mathcal {R}}}_1 - \frac{1}{2}{{\mathcal {R}}}_T{{\mathcal {R}}}_VX^\textrm{T}Z \right) \left( I+\frac{t}{2}{{\mathcal {R}}}_VX^\textrm{T}Z\right) X^\textrm{T}X \end{aligned}$$

and ${{\mathcal {R}}}_T^\textrm{Full} = {{\mathcal {R}}}_T X^\textrm{T}X$, we have that

$$\begin{aligned} \left\{ \begin{array}{l} \displaystyle W_1(t) = \left( I + \frac{t}{2}{{\mathcal {R}}}_V X^\textrm{T}Z\right) X^\textrm{T}X \left( I + \frac{t}{2}W_2(t)\right) + O(\Vert {{\mathcal {R}}}_V\Vert ^2), \quad \text{ and } \\ \displaystyle W_2(t) = \left( I - \frac{t}{2}\left( -Z^\textrm{T}X - \frac{t}{2}{{\mathcal {R}}}_1^\textrm{Full} + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2}\right) \right) ^{-1} \left( -Z^\textrm{T}X - \frac{t}{2}{{\mathcal {R}}}_1^\textrm{Full} + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2}\right) \end{array} \right. \end{aligned}$$

(A1)

Proof

From Theorem 1, we have that $X^{+}(t) = X + tUW(t)$, with W(t) solution of

$$\begin{aligned} \left( I-\frac{t}{2}V^\textrm{T}U\right) W(t) = V^\textrm{T}X. \end{aligned}$$

(A2)

Now, using definition of ${{\mathcal {R}}}_V$ and ${{\mathcal {R}}}_T$, we have, from definition of ${{\mathcal {R}}}_1^\textrm{Full}$, that

$$\begin{aligned} \small V^\textrm{T} U = \left( \begin{array}{cc} (I-XX)X^\textrm{T}Z &{} X^\textrm{T}X \\ \displaystyle (-Z+\frac{1}{2}X{{\mathcal {R}}}_T)^\textrm{T}(Z-XX^\textrm{T}Z) &{}\displaystyle -Z^\textrm{T}X+\frac{{{\mathcal {R}}}_T^\textrm{Full}}{2} \end{array}\right) = \left( \begin{array}{cc} {{\mathcal {R}}}_VX^\textrm{T}Z &{} X^\textrm{T}X \\ \displaystyle -{{\mathcal {R}}}_1^\textrm{Full} &{} \displaystyle -Z^\textrm{T}X+\frac{{{\mathcal {R}}}_T^\textrm{Full}}{2} \end{array}\right) \end{aligned}$$

and

$$\begin{aligned} V^\textrm{T} X = \left( \begin{array}{cc} X^\textrm{T}X \\ \displaystyle -Z^\textrm{T}X + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2} \end{array}\right) . \end{aligned}$$

Thus, linear system (A2) leads to

$$\begin{aligned} \left\{ \begin{array}{l} \displaystyle W_1(t) = \left( I - \frac{t}{2}{{\mathcal {R}}}_V X^\textrm{T}Z\right) ^{-1} X^\textrm{T}X\left( I + \frac{t}{2}W_2(t)\right) , \quad \text{ and } \\ \displaystyle \frac{t}{2} {{\mathcal {R}}}_1^\textrm{Full} W_1(t) + \left( I - \frac{t}{2}\left( -Z^\textrm{T}X + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2}\right) \right) W_2(t) = -Z^\textrm{T}X + \frac{1}{2} {{\mathcal {R}}}_T^\textrm{Full} \end{array} \right. \end{aligned}$$

Now, by assuming that $\Vert \frac{t}{2}{{\mathcal {R}}}_V X^\textrm{T}Z\Vert < 1$, from Bannach’s Lemma (Golub and Van Loan 1996), $(I - \frac{t}{2}{{\mathcal {R}}}_V X^\textrm{T}Z)^{-1} = I + \frac{t}{2}{{\mathcal {R}}}_V X^\textrm{T}Z + O(\Vert {{\mathcal {R}}}_V\Vert ^2)$. Then, previous linear system becomes

$$\begin{aligned} \left\{ \begin{array}{l} \displaystyle W_1(t) = \left( I + \frac{t}{2}{{\mathcal {R}}}_VX^\textrm{T}Z\right) X^\textrm{T}X\left( I + \frac{t}{2}W_2(t)\right) + O(\Vert {{\mathcal {R}}}_V\Vert ^2), \\ \displaystyle W_2(t) = \left( I - \frac{t}{2}\left( -Z^\textrm{T}X - \frac{t}{2}{{\mathcal {R}}}_1^\textrm{Full} + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2}\right) \right) ^{-1}\left( -Z^\textrm{T}X - \frac{t}{2}{{\mathcal {R}}}_1^\textrm{Full} + \frac{{{\mathcal {R}}}_T^\textrm{Full}}{2}\right) \end{array} \right. \end{aligned}$$

from where the proof follows. $\square$

It is worth observing that when ${{\mathcal {R}}}_V =0$, Theorem 3 is reduced to Theorem 1. In addition, since ${{\mathcal {R}}}_V \approx 0$ (since X is practically feasible), the term $O(\Vert {{\mathcal {R}}}_V\Vert ^2)$ can be disregarded in (A1).