Approximation of Matrix Functions Arising in Physics and Network Science: Theoretical and Computational Aspects

Many applications in physics and network science require the computation of quantities related to certain matrix functions. In many cases, a straightforward way to proceed is by diagonalization. However, the cost of this method becomes prohibitive for large dimensions. We aim to provide techniques with much lower computational cost that exploit the (approximate) sparse structure of the matrices involved. Spectral projectors associated with Hermitian matrices play a key role in applications such as electronic structure computations. Linear scaling methods in the gapped case are based on the fact that the entries of such projectors decay rapidly with respect to some sparsity patterns. The relation with the sign function, together with an integral representation of the latter, is used to obtain asymptotically optimal decay bounds. The influence of isolated extremal eigenvalues on the decay properties is also investigated. Using similar techniques, we extend our results to related matrix functions. Another problem of interest is the computation of the trace of a matrix function. In certain situations, combining probing methods based on graph colorings with stochastic trace estimation techniques can yield accurate approximations at moderate cost. Until recently, such methods had not been thoroughly analyzed, however, but were rather used as efficient heuristics by practitioners. We perform a detailed analysis of stochastic probing methods and, in particular, expose conditions under which the expected approximation error in the stochastic probing method scales more favorably with the dimension of the matrix than the error in non-stochastic probing. A quantity that often appears in quantum physics and network science is the von Neumann entropy, expressed as the trace of a matrix function. Probing techniques or stochastic trace estimators can be used to obtain approximations of the entropy. The computation of the several quadratic forms and matrix-vector products involved can be carried efficiently by using polynomial and rational Krylov subspace methods. For the probing approach, theoretical bounds and heuristic estimates are provided for the error on the entropy, which can be used to select the number of quadratic forms required to reach a certain accuracy. Moreover, a posteriori error bounds are given for the Krylov approximations. Our results are validated by several numerical experiments on a number of test problems arising in network analysis.