Fast way of computing the action of a matrix power on a vector

Question

For integer $k>0$, it is well-known that one can use binary exponentiation to evaluate the matrix power $\mathbf A^k$, where $\mathbf A$ is an $n\times n$ matrix.

However, it is not clear to me if that method can be readily modified to evaluate $\mathbf A^k \mathbf v$, where $\mathbf v$ is a vector of dimension $n$. Of course, one could just form the power as usual, and then apply it to the vector, but in cases like this (e.g. if the matrix is very large), one would wish to avoid forming a whole matrix. (At worst, one could use a few vectors of comparable size.)

Is there a good way to compute the action of a matrix power on a vector for integer exponents? My attempts to search the literature have come up empty, but perhaps I am not using the right keywords.

Solutions, or at least pointers to the literature, would be very much appreciated.

Answer 1

It's going to be hard to beat the two naive approaches, iterated squaring and repeated application associating from the right $A(A(A(\cdots Av)))$. Which one of the two is better depends on the value of $k$, sparsity, etc. as noted in the comments.

If you are fine with approximations, there might be other options; for instance, you could:

(1) treat your problem as a general $f(A)v$, where $f$ is a scalar function. This opens up a few strategies, for instance the diagonalization or Schur--Parlett method: a change of basis to reduce to diagonal or triangular form + solving explicitly the triangular/diagonal problem. This might be convenient for dense $A$ and huge $k$ since its complexity does not depend on $k$ (apart from computing $k$th powers of scalars, which can usually be done in constant time). Or the Arnoldi method for sparse $A$, which similarly has a very mild dependence on the exponent $k$. These are treated on Higham's book Functions of matrices.

(2) approximate $A$ with a rank-$r$ matrix $UV^T$, for instance with randomized methods, and then use $(UV^T)^{k} = U(V^TU)^{k-1}V^T$, which leaves you with powers of a $r\times r$ matrix rather than a $n\times n$ one. The techniques for this random approximation are sketched (pun intended) for instance in the popular review article https://arxiv.org/abs/2002.01387 . Note that the Arnoldi method suggested above is also essentially a rank-$r$ approximation computed in a peculiar way.