The semi-explicit nonsmooth Newmark-β time integrator for robust unilateral contact in dynamic fragmentation simulations

Numerical simulations of solids undergoing dynamic fragmentation, a problem characterized by dynamic fracture and dense contacts, require accurately capturing the transition from a solid continuum to a collection of interacting fragments. We use the finite-element method with the extrinsic cohesive zone model for fracture. For contact, conventional penalty-based methods often exhibit numerical instabilities in dynamic collision-rich settings. To address this, we adapt and validate a novel semi-explicit time-integration scheme: the Nonsmooth Newmark-β (NSN) method for unilateral contact. Based on the Nonsmooth Contact Dynamics (NSCD) method, this formulation strongly enforces contact constraints at the velocity level. Within this scheme, the bulk dynamics are non-impulsive and integrated explicitly with second-order accuracy, and the fracture model allows displacement discontinuities and integrates the cohesive softening explicitly. Contact, treated rigorously as nonsmooth, is integrated implicitly with first-order accuracy. Benchmark tests demonstrate that the NSN scheme achieves accuracy comparable to established nonsmooth methods, such as the semi-explicit CD-Lagrange and implicit Moreau–Jean schemes. Moreover, it outperforms penalty-based approaches by orders of magnitude. Although the NSN method incurs a higher per-step computational cost, its enhanced stability allows for significantly larger time steps. Consequently, for 1D benchmarks, overall computational efficiency is comparable to or better than that of purely explicit approaches. We applied this framework to 1D fragmentation under free and confined expansion. Results reveal that confinement shifts the fracture energy budget from local fragment kinetic energy to larger-scale global system kinetic energy. Additionally, we found, counterintuitively, that, compared to fully elastic contact, adding contact dissipation reduces fracture energy yet increases the final fragment count. It occurs because such dissipation reduces the vibration within damaged fragments, allowing cleaner stress-wave propagation and better damage localization, driving cracks to full separation rather than distributing damage. These results establish the NSN scheme as a robust tool for generating high-fidelity fragmentation statistics.