Solving sparse linear systems from discretized PDEs is challenging. Direct solvers have in many cases quadratic complexity (depending on geometry), while iterative solvers require problem dependent preconditioners to be robust and efficient. Approximate factorization preconditioners, such as incomplete LU factorization, provide cheap approximations to the system matrix. However, even a highly accurate preconditioner may have deteriorating performance when the condition number of the system matrix increases. By increasing the accuracy on low-frequency errors, we propose a novel hierarchical solver with improved robustness with respect to the condition number of the linear system. This solver retains the linear computational cost and memory footprint of the original algorithm.