Numerical methods for exascale magnetohydrodynamics simulations in general relativity