Convergence Analysis of Markov Chain Monte Carlo Linear Solvers Using Ulam--von Neumann Algorithm

Abstract

The convergence of Markov chain-based Monte Carlo linear solvers using the Ulamvon Neumann algorithm for a linear system of the form x = Hx + b is investigated in this paper. We analyze the convergence of the Monte Carlo solver based on the original Ulam-von Neumann algorithm under the conditions that H < 1 as well as ρ(H) < 1, where ρ(H) is the spectral radius of H. We find that although the Monte Carlo solver is based on sampling the Neumann series, the convergence of Neumann series is not a sufficient condition for the convergence of the Monte Carlo solver. Actually, properties of H are not the only factors determining the convergence of the Monte Carlo solver; the underlying transition probability matrix plays an important role. An improper selection of the transition matrix may result in divergence even though the condition H < 1 holds. However, if the condition H < 1 is satisfied, we show that there always exist certain transition matrices that guarantee convergence of the Monte Carlo solver. On the other hand, if ρ(H) < 1 but H ≥ 1, the Monte Carlo linear solver may or may not converge. In particular, if the row sum N j=1 |H ij | > 1 for every row in H or, more generally, ρ(H +) > 1, where H + is the nonnegative matrix where H + ij = |H ij |, we show that transition matrices leading to convergence of the Monte Carlo solver do not exist. Finally, given H and a transition matrix P , denoting the matrix H * via H * ij = H 2 ij /P ij , we find that ρ(H *) < 1 is a necessary and sufficient condition for convergence of the Markov chain-based Monte Carlo linear solvers using the Ulam-von Neumann algorithm.

Multiway Monte Carlo Method for Linear Systems

“…To tackle this problem, we propose a multi-way Markov random walk which uses multiple transition matrices. At each step of the random walk, the transition matrix is constructed in a way akin to the Monte Carlo Almost Optimal (MAO) framework [8,12]. We prove that under this type of random walk, the new method always converges when ρ(H + ) < 1, where H + is the nonnegative matrix as H + ij = |H ij |.…”

“…Examples of these works are Sequential Monte Carlo method [11] and synthetic-acceleration method [10]. Also there are a variety of studies of the parallel implementation [6,14,1], real world application [17,2], convergence analysis [12], and spectral analysis [15].…”

A highly parallel algorithm for computing the action of a matrix exponential on a vector based on a multilevel Monte Carlo method

“…The method can be understood formally as a procedure consisting in a Monte Carlo sampling of the Neumann series of the inverse of the matrix. The convergence of the method was rigorously established in [29], and improved further more recently (see for instance [12], and [19] just to cite a few references).…”

A Monte Carlo method for computing the action of a matrix exponential on a vector

“…Mathematically, the method can be seen in a way as a Monte Carlo sampling of the Neumann series of the matrix. The convergence of the method was rigorously established in [26], and improved further more recently (see for instance [15], and [7] just to cite a few references). More recently, and for the specific case of computing the action of a Hermitian matrix exponential over a vector, which is of interest in Quantum Physics, it has been proposed in [32] an efficient algorithm based on a novel randomized linear algebra technique known in the literature as the Nyström method.…”