  A new version of a numerical model of stellar differential rotation based on
mean-field hydrodynamics is presented and tested by computing the differential
rotation of the Sun. The model is then applied to four individual stars
including two moderate and two fast rotators to reproduce their observed
differential rotation quite closely. A series of models for rapidly rotating
($P_{rot} = 1$ day) stars of different masses and compositions is generated.
The effective temperature is found convenient to parameterize the differential
rotation: variations with metallicity, that are quite pronounced when the
differential rotation is considered as a function of the stellar mass, almost
disappear in the dependence of differential rotation on temperature. The
differential rotation increases steadily with surface temperature to exceed the
largest differential rotation observed to date for the hottest F-stars we
considered. This strong differential rotation is, however, found not to be
efficient for dynamos when the efficiency is estimated with the standard
$C_\Omega$-parameter of dynamo models. On the contrary, the small differential
rotation of M-stars is the most dynamo-efficient. The meridional flow near the
bottom of the convection zone is not small compared to the flow at the top in
all our computations. The flow is distributed over the entire convection zone
in slow rotators but retreats to the convection zone boundaries with increasing
rotation rate, to consist of two near-boundary jets in rapid rotators. The
implications of the change of the flow structure for stellar dynamos are
briefly discussed.
