A new theoretical approach for turbulent flows based on Lie-group analysis is presented. It unifies a large set of ‘solutions’ for the mean velocity of stationary parallel turbulent shear flows. These results are not solutions in the classical sense but instead are defined by the maximum number of possible symmetries, only restricted by the flow geometry and other external constraints. The approach is derived from the Reynolds-averaged Navier–Stokes equations, the fluctuation equations, and the velocity product equations, which are the dyad product of the velocity fluctuations with the equations for the velocity fluctuations. The results include the logarithmic law of the wall, an algebraic law, the viscous sublayer, the linear region in the centre of a Couette flow and in the centre of a rotating channel flow, and a new exponential mean velocity profile not previously reported that is found in the mid-wake region of high Reynolds number flat-plate boundary layers. The algebraic scaling law is confirmed in both the centre and the near-wall regions in both experimental and DNS data of turbulent channel flows. In the case of the logarithmic law of the wall, the scaling with the distance from the wall arises as a result of the analysis and has not been assumed in the derivation. All solutions are consistent with the similarity of the velocity product equations to arbitrary order. A method to derive the mean velocity profiles directly from the two-point correlation equations is shown.