The development of fast numerical methods for multilevel radiative transfer (RT) applications often leads to important breakthroughs in astrophysics, becausethey allow the investigation of problems that could not be properly tackled using the methods previously available. Probably, themost familiar example is the so-called Multilevel Accelerated Λ-Iteration (MALI) technique of Rybicki & Hummerfor the case of a local approximate operator,which is based on Jacobi iteration. However, there are superior operator-splitting methods, based on Gauss-Seidel (GS) and Successive Overrelaxation (SOR) iteration, which provide a dramatic increase in the speed with which non-LTE multilevel transfer problems can be solved in one, two and three-dimensional geometries. Such RT methods, which were introduced by Trujillo Bueno & Fabiani Bendicho ten years ago, are the main subject ofthe first part of this paper. We show in some detail how they can be applied for solving multilevel RT problems in spherical geometry, for both atomic and molecular line transitions. The second part of the article addresses the issue of the calculation of the molecular number densities when the approximation of instantaneous chemical equilibrium turns out to be inadequate, which happens to be the case wheneverthe dynamical time scales of the astrophysical plasma under considerationare much shorter than the time needed by the molecules to form.