CFD simulations of gas-liquid flows have been carried out in a full 3D, unsteady, Euler/Euler framework by means of the commercial software FLUENT(R) and discussed in connection with laboratory scale bubble columns, with a particular interest in churn-turbulent flows. The importance of the choice of third order numerical schemes for such simulations was illustrated, as significant errors due to numerical diffusion can be induced by other less-stringent numerical set-ups. The available choices for turbulence modeling in RANS approach were tested following 9 different options, i.e., three possible formulations of the k-epsilon model (Standard, RNG, Realizable) combined with three different modalities to account for gas-phase effects (Dispersed, Dispersed+Bubble Induced Turbulence, Per-Phase). The Standard and Realizable versions of the k-epsilon model were unable to reproduce the expected gulf-stream patterns of bubble columns. However, the RNG k-epsilon model has exhibited much better descriptions of the flow features. The inclusion of bubble-induced terms in the turbulence equations led to minor impact on the performance of RNG k-epsilon models. The superior performance of RNG k-epsilon models essentially ties in the better depiction of the turbulent dissipation rate and turbulent viscosity. Thus, the application of RNG k-epsilon models would be also recommended for a further implementation of bubble population balance models. (C) 2009 Elsevier Ltd. All rights reserved.