The specific heat capacity ( \({C}_{p}\) ) of submerged arc welding fluxes plays a critical role in determining heat-storage capacity and, in conjunction with thermal conductivity and density, collectively governs heat transfer kinetics during welding. However, the microscopic factors of \({C}_{p}\) variation in multicomponent silicate fluxes remain unclear. In this study, molecular dynamics (MD) simulations combined with vibrational density of states (VDOS) analysis and quantum heat-capacity correction have been employed to investigate the effect of ZrO2 on \({C}_{p}\) of CaF2–SiO2–CaO–ZrO2 slags. The results confirm a reduction in \({C}_{p}\) with increasing ZrO2 content. Structural analysis reveals that replacement of Si–O–Si by more rigid Zr–O–Si bonds imposes local geometric constraints on the silicate network. This local structural tightening drives Si-related vibrational modes toward higher frequencies and the local cancellation of oxygen vibrational activity. This vibrational reorganization inherently suppresses the intrinsic vibrational heat capacity ( \({C}_\text{vib}\) ), thereby serving as a rigorous and quantifiable microscopic fingerprint that governs the \({C}_{p}\) degradation. This work establishes a correlation between atomic-scale vibrational characteristics and \({C}_{p}\) , advancing the predictive capability of thermophysical behaviors in sophisticated welding fluxes.