Aviscoelastic model (a mixing-length model with relaxation) is developed to investigate the effect of turbulent advection on the mean flow perturbation and the drag force induced by turbulent shear flow over an undulating surface. The relaxation term is proportional to the ratio of eddy turnover time to travelling time; accordingly, near the surface, the relaxation model reduces to an eddy-viscosity or mixing-length model, while far from the surface it reduces to a rapid-distortion model.The linear governing equations are transformed into streamline coordinates and solved through matched asymptotic expansions. According to order-of-magnitude estimates in Belcher, Newley & Hunt (1993), the drag force contributed by nonlinear shear stress is of the same order as that contributed by asymmetric pressure arising from the leeward thickening of the perturbed boundary layer. The nonlinear analysis in the present model confirms this estimate in most cases. Our analytical results show a dip in shear stress at the interface between the inner and outer layers and provide evidence that this dynamical feature is related to eddy advection. Numerical calculation using a shooting method gives results that compare well with the analysis.