The discrete-element method was used to simulate drained triaxial compression of large-scale, polydisperse numerical samples at a range of void ratios while tracing all relevant energy components. The frictional dissipation and boundary work are almost equal regardless of sample density. The volumetric work reaches a steady value at large strain. However, the distortional work increases continually as sample deformation continues post-critical state. There is a preferential orientation for frictional dissipation at around 45° to the major principal stress direction. This matches the orientation at which there is the largest number of sliding contacts. The work equations, which are fundamental in most commonly used constitutive models, are linear when plotted against deviatoric strain. The modified Cam Clay work equation substantially over-predicts the frictional dissipation for dense samples. An alternative, thermodynamically consistent work equation gives a much better description of frictional dissipation and is therefore recommended to ensure accuracy in modelling.