We report a robust numerical method that allows the simulation of fibre suspensions, together with shear-induced fibre migration, to be carried out at a large volume fraction of the fibres. The method is a combination of the Brownian configuration field (BCF) technique and the adaptive viscosity split stress (AVSS) finite element formulation. The fibre suspension model is adapted from Folgar and Tucker [E.P.. Folgar, C.L. Tucker, J. Reinforced Plastics and Composites 3 (1984) 98-119] (for the structure evolution) and Phan-Thien and Graham [N. Phan-Thien, A.L. Graham, Rheol. Acta 30 (1991) 44-57] (for the stress rule), together with a modified version of the shear-induced migration model of Phillips et al. [R.J. Phillips, R.C. Armstrong, R.A. Brown, A.L. Graham, J.R. Abott, Phys. Fluids A4 (1992) 30-40]. The implementation was tested in the circular Couette and the plane Poiseuille flow. The numerical data compared well with the experimental results of Mondy et al. [L.A. Mondy, H. Brenner, S.A. Altobelli, J.R. Abbott, A.L. Graham, J. Rheol. 38 (1994) 444-452]. (C) 2000 Elsevier Science B.V. All rights reserved.