We present an efficient computational method for calculating the binding energies of the bound states of ultracold alkali-metal dimers in the presence of magnetic fields. The method is based on propagation of coupled differential equations and does not use a basis set for the interatomic distance coordinate. It is much more efficient than the previous method based on a radial basis set and allows many more spin channels to be included. This is particularly important in the vicinity of avoided crossings between bound states. We characterize a number of different avoided crossings in Cs-2 and compare our converged calculations with experimental results. Small but significant discrepancies are observed in both crossing strengths and level positions, especially for levels with l symmetry (rotational angular momentum L=8). The discrepancies should allow the development of improved potential models in the future.