Abstract Trace impurity transport is studied with the flux-driven gyrokinetic GYSELA code [V. Grandgirard et al., emphComp. Phys. Commun. bf 207, 35 (2016)]. A reduced and linearized multi-species collision operator has been recently implemented, so that both neoclassical and turbulent transport channels can be treated self-consistently on an equal footing. In the Pfirsch-Schlüter regime likely relevant for tungsten, the standard expression of the neoclassical impurity flux is shown to be recovered from gyrokinetics with the employed collision operator. Purely neoclassical simulations of deuterium plasma with trace impurities of helium, carbon and tungsten lead to impurity diffusion coefficients, inward pinch velocities due to density peaking, and thermo-diffusion terms which quantitatively agree with neoclassical predictions and NEO simulations [E. Belli et al., emphPlasma Phys. Control. Fusion textbf54, 015015 (2012)]. The thermal screening factor appears to be less than predicted analytically in the Pfirsch-Schlẗer regime, which can be detrimental to fusion performance. Finally, self-consistent nonlinear simulations have revealed that the tungsten impurity flux is not the sum of turbulent and neoclassical fluxes computed separately, as usually assumed. The synergy partly results from the turbulence-driven in-out poloidal asymmetry of tungsten density. This result puts forward the need for self-consistent simulations of impurity transport, i.e. including both turbulence and neoclassical physics, in view of quantitative predictions for ITER.