Hyperspectral super-resolution (HSR) is a problem of recovering a high-spectral-spatial-resolution image from a multispectral measurement and a hyperspectral measurement, which have low spectral and spatial resolutions, respectively. We consider a low-rank structured matrix factorization formulation for HSR, which is a non-convex large-scale optimization problem. Our contributions contain both computational and theoretical aspects. On the computational side, we develop three inexact block coordinate descent (BCD) schemes that are empirically found to run many times faster than a state-of-the-art method, which uses exact BCD. We achieve this by applying concepts in the proximal gradient (PG) and Frank-Wolfe (FW) methods and by exploiting the HSR problem structures. On the theoretical side, we show that these inexact BCD schemes guarantee convergence to a stationary point. In particular, the convergence result for a hybrid PG-FW inexact BCD scheme is new.