This talk is concerned with representation issues associated with the numerical solution of a unified mathematical model of continuum mechanics, due to Godunov, Peshkov, and Romenski, which can describe ideal fluids, viscous fluids and elastoplastic solids as special cases of a general continuum. The different regimes are characterized solely by the choice of material parameters and the resulting PDE system is of hyperbolic nature, with clearly defined finite wave speeds, in contrast to the standard formulation of viscous fluids via the Navier--Stokes equations. The description of such a general continuum hinges on the evolution of a matrix-valued field called distortion, which is a generalization of the inverse deformation gradient in solid mechanics. In the fluid regime, this quantity can no longer be recovered as a gradient of displacements and encodes very rich information, in particular due to the different orientations that ideal fluid parcels can be found in. The fine features of the distortion field can be challenging (or outright impossible) to resolve with standard well-tested Finite Volume methods. Degenerate situations are routinely encountered where unphysical states are generated simply as a result of taking a convex combination of two data points. We show how changing to an alternative representation of the same distortion field, obtained via polar decomposition, can be used to solve such discretization issues. Instead of the original PDE system, one can instead evolve the rotational and stretch components of the distortion matrix separately, which allows the description of the rotational components through a quaternion-valued partial differential equation. We discuss the peculiarities of quaternion PDEs and some of the discretization strategies that they enable. We present numerical examples of high-Reynolds number simulations which could not be carried out with the previous formulation of the model.