The diffusive strip method (DSM) is a near-exact numerical method for mixing computations initially developed in two dimensions (Meunier & Villermaux, J. Fluid Mech., vol. 662, 2010, pp. 134–172). The method, which consists of following stretched material lines to compute the resulting scalar field a posteriori, is extended here to three-dimensional flows. We describe the procedure and its three-dimensional peculiarity, which relies on the Lagrangian advection of a triangulated surface from which the stretching rate is extracted to infer the scalar field. The method is first validated at moderate Péclet number against a classical pseudospectral method solving the advection–diffusion equation for a Batchelor vortex, and then applied to a simple Taylor–Couette experimental configuration with non-rotating boundary conditions at the top-end disk, bottom-end disk and outer cylinder. This motion, producing an elaborate although controlled steady three-dimensional flow, relies on Ekman pumping arising from the rotation of the inner cylinder. A recurrent two-cell structure is separated by the horizontal mid-plane and formed by stream tubes shaped as nested tori under laminar flow conditions. A scalar blob in the flow experiences a Lagrangian oscillating dynamics undergoing stretchings and compressions, driving the mixing process. The DSM enables the calculation of the blob elongation and scalar concentration distributions through a single variable computation along the advected blob surface, capturing the rich evolution observed in the experiments. Interestingly, the mixing process in this axisymmetric and steady three-dimensional flow leads to a linear growth of surfaces in time similar to the one obtained in a two-dimensional shear. The potentialities, limits and extension of the method to more general flows are finally discussed.