  In this paper we investigate the effectiveness of direct statistical
simulation (DSS) for two low-order models of dynamo action. The first model,
which is a simple model of solar and stellar dynamo action, is third-order and
has cubic nonlinearities whilst the second has only quadratic nonlinearities
and describes the interaction of convection and an aperiodically reversing
magnetic field. We show how DSS can be utilised to solve for the statistics of
these systems of equations both in the presence and the absence of stochastic
terms, by truncating the cumulant hierarchy at either second or third order. We
compare two different techniques for solving for the statistics, timestepping
-- which is able to locate only stable solutions of the equations for the
statistics and direct detection of the fixed points. We develop a complete
methodology and symbolic package in Python for deriving the statistical
equations governing the Low-order dynamic systems in cumulant expansions. We
demonstrate that although direct detection of the fixed points is efficient and
accurate for DSS truncated at second order, the addition of higher order terms
leads to the inclusion of many unstable fixed points that may be found by
direct detection of the fixed point by iterative methods. In those cases
timestepping is a more robust protocol for finding meaninful solutions to DSS.
