Multi-variable cell population balance models represent the most accurate and general way of describing the complicated phenomena associated with cell growth, substrate consumption and product formation due to the level of detail included in them. Therefore, the ability to solve and understand such models is of fundamental importance in predicting and/or controlling cell growth in processes of biotechnological interest. However, due to the fact that such models typically consist of first-order, partial integro-differential equations coupled in a nonlinear fashion with ordinary integro-differential equations, their solution poses a serious challenge. In this work, we have developed several finite difference algorithms for the solution of the problem in its most general formulation (i.e. for any set of single-cell physiological state functions). The validity of the developed algorithms was verified by comparing their results with those of three specific test problems for which several solution characteristics are known. Moreover, the numerical schemes were compared to each other with respect to their key numerical features, such as stability, accuracy and computational speed. Solutions of the cell population balance model with up to three state variables were obtained using a Pentium II 450 MHz PC in tractable CPU times.