Prethermal gauge structure and surface growth in mathbb{Z}₂ lattice gauge theories

Universal aspects of thermalization in interacting many-body systems are challenging to derive microscopically, especially in kinetically constrained models, yet their numerical study beyond (1+1)D remains notoriously difficult. Here, we numerically study the mean-field dynamics of a (2+1)D spin system with thousands of spins and show that experimentally-feasible two-body Ising interactions can stabilize a prethermal mathbb{Z}₂ lattice gauge structure with dynamical matter, manifested by a separation of timescales with a stable gauge-invariant plateau. Eventually, the metastable prethermal mathbb{Z}₂ gauge structure breaks down via a proliferation of Gauss' law defects, similar to bubble formation in false vacuum decay. In this regime, we discover spatio-temporal correlations described by a non-linear surface growth consistent with the (1+1)D Kardar-Parisi-Zhang (KPZ) universality class, revealing a previously hidden feature in the thermalization of multi-point correlators. We benchmark our results in small systems against semi-classical discrete time Wigner approximation (DTWA) and exact diagonalization (ED), where the breakdown of DTWA signals the emergence of an extensive number of local symmetries that strongly influence the thermalization pathway. Our model provides a testbed for quantum simulators and is directly implementable in large-scale arrays of Rydberg atoms.