Summary: In many engineering applications it is required to compute the dominant subspace of a matrix \(A\) of dimension \(m\times n\), with \(m\gg n\). Often the matrix \(A\) is produced incrementally, so all the columns are not available simultaneously. This problem arises, e.g., in image processing, where each column of the matrix \(A\) represents an image of a given sequence leading to a singular value decomposition-based compression [S. Chandrasekaran, B.S. Manjunath, Y.F. Wang, J. Winkeler and H. Zhang, An eigenspace update algorithm for image analysis, Graphical Models and Image Process. 59 (5) 321–332 (1997)]. Furthermore, the so-called proper orthogonal decomposition approximation uses the left dominant subspace of a matrix \(A\) where a column consists of a time instance of the solution of an evolution equation, e.g., the flow field from a fluid dynamics simulation. Since these flow fields tend to be very large, only a small number can be stored efficiently during the simulation, and therefore an incremental approach is useful [P. Van Dooren, in: Griffiths, D. F. (ed.) et al., Numerical analysis 1999. Proceedings of the 18th Dundee biennial conference, Univ. of Dundee, GB, June 29th–July 2nd, 1999. Boca Raton, FL: Chapman & Hall/CRC. Chapman Hall/CRC Res. Notes Math. 420, 231–247 (2000)].
In this paper, an algorithm for computing an approximation of the left dominant subspace of size \(k\) of \(A\in \mathbb R^{m\times n}\), with \(k\ll m,n\), is proposed requiring at each iteration O\((mk+k^{2})\) floating point operations. Moreover, the proposed algorithm exhibits a lot of parallelism that can be exploited for a suitable implementation on a parallel computer.