  Nonlinear mean-field models of the solar dynamo show long-term variability,
which may be relevant to different states of activity inferred from long-term
radiocarbon data. This paper is aimed to probe the dynamo hysteresis predicted
by the recent mean-field models of Kitchatinov \& Olemskoy (2010) with direct
numerical simulations. We perform three-dimensional simulations of large-scale
dynamos in a shearing box with helically forced turbulence. As initial
condition, we either take a weak random magnetic field or we start from a
snapshot of an earlier simulation. Two quasi-stable states are found to coexist
in a certain range of parameters close to the onset of the large-scale dynamo.
The simulations converge to one of these states depending on the initial
conditions. When either the fractional helicity or the magnetic Prandtl number
is increased between successive runs above the critical value for onset of the
dynamo, the field strength jumps to a finite value. However, when the
fractional helicity or the magnetic Prandtl number is then decreased again, the
field strength stays at a similar value (strong field branch) even below the
original onset. We also observe intermittent decaying phases away from the
strong field branch close to the point where large-scale dynamo action is just
possible. The dynamo hysteresis seen previously in mean-field models is thus
reproduced by 3D simulations. Its possible relation to distinct modes of solar
activity such as grand minima is discussed.
