Two-dimensional axisymmetric Eulerian/Eulerian simulations of two-phase (gas/liquid) transient flow were performed using a multiphase flow algorithm based on the finite-volume method. These numerical simulations cover laboratory scale bubble columns of different diameters, operated over a range of superficial gas velocities ranging from the bubbly to the churn turbulent regime. The bubble population balance equation (BPBE) is implemented in the two-fluid model that accounts for the drag force and employs the modified k-epsilon turbulence model in the liquid phase. Several available bubble breakup and coalescence closures are tested. Quantitative agreements between the experimental data and simulations are obtained for the time-averaged axial liquid velocity profiles, as well as for the kinetic energy profiles, only when model predicted breakup rate is increased by a factor of ten to match the coalescence rate. The calculated time-averaged gas holdup profiles deviate in shape from the measured ones and suggest that full three-dimensional simulation is needed. Implementation of BPBE leads to better agreement with data, especially in the chum-turbulent flow regime, compared to the simulation based on an estimated constant mean bubble diameter. Differences in the predicted interfacial area density, with and without BPBE implementation, are significant. The choice of bubble breakup and coalescence closure does not have a significant impact on the simulated results as long as the magnitude of breakup is increased tenfold. (C) 2004 Elsevier Ltd. All rights reserved.