  Increasing evidence is becoming available about not only the surface
differential rotation of rapidly rotating cool stars but, in a small number of
cases, also about temporal variations, which possibly are analogous to the
solar torsional oscillations. Given the present difficulties in resolving the
precise nature of such variations, due to both the short length and poor
resolution of the available data, theoretical input is vital to help assess the
modes of behaviour that might be expected, and will facilitate interpretation
of the observations. Here we take a first step in this direction by studying
the variations in the convection zones of such stars, using a two dimensional
axisymmetric mean field dynamo model operating in a spherical shell in which
the only nonlinearity is the action of the azimuthal component of the Lorentz
force of the dynamo generated magnetic field on the stellar angular velocity.
We consider three families of models with different depths of dynamo-active
regions. For moderately supercritical dynamo numbers we find torsional
oscillations that penetrate all the way down to the bottom of the convection
zones, similar to the case of the Sun. For larger dynamo numbers we find
fragmentation in some cases and sometimes there are other dynamical modes of
behaviour, including quasi-periodicity and chaos. We find that the largest
deviations in the angular velocity distribution caused by the Lorentz force are
of the order of few percent, implying that the original assumed `background'
rotation field is not strongly distorted.
