Numerical integration methods for stochastic volatility models in financial markets are discussed. We concentrate on two classes of stochastic volatility models where the volatility is either directly given by a mean-reverting CEV process or as a transformed Ornstein-Uhlenbeck process. For the latter, we introduce a new model based on a simple hyperbolic transformation. Various numerical methods for integrating mean-reverting CEV processes are analysed and compared with respect to positivity preservation and efficiency. Moreover, we develop a simple and robust integration scheme for the two-dimensional system using the strong convergence behaviour as an indicator for the approximation quality. This method, which we refer to as the IJK (137) scheme, is applicable to all types of stochastic volatility models and can be employed as a drop-in replacement for the standard log-Euler procedure.