  Calculating cost-effective solutions to particle dynamics in viscous flows is
an important problem in many areas of industry and nature. We implement a
second-order symmetric splitting method on the governing equations for a rigid
spheroidal particle model with torques, drag and gravity. The method splits the
operators into a vector field that is conservative and one that takes into
account the forces of the fluid. Error analysis and numerical tests are
performed on perturbed and stiff particle-fluid systems. For the perturbed
case, the splitting method greatly improves the solution accuracy, when
compared to a conventional multi-step method, and the global error behaves as
$\mathcal{O}(\varepsilon h^2)$ for roughly equal computational cost. For stiff
systems, we show that the splitting method retains stability in regimes where
conventional methods blow up. In addition, we show through numerical
experiments that the global order is reduced from
$\mathcal{O}(h^2/\varepsilon)$ in the non-stiff regime to $\mathcal{O}(h)$ in
the stiff regime.
