Quasi-brittle materials create intricate crack patterns on their surface and internal three-dimensional crack networks when exposed to thermal shock conditions. Their cracking behavior is notably impacted by their quenching temperature and size. In this paper, a coupled three-dimensional phase field-cohesive zone model that accounts for thermoelastic fracture is developed to reproduce the quenching and cracking experiments of ceramic balls. The accuracy of the phase field model is demonstrated by simulating the phase field of a single-element cooling process and comparing it with the analytical solution. The quench fracture process of the 3D ceramic ball is then simulated. The phase field simulations show that the surface crack morphology during the quenching process is somewhat random. Still, the crack distribution density is determined by the quenching temperature difference and ceramic ball size. Furthermore, a threshold size is present during the quenching process. When the ceramic ball's radius falls below this threshold, the surface crack pattern no longer appears during quenching. Consequently, the phase-field-cohesive zone model represents an effective means of predicting three-dimensional fracture processes in quasi-brittle materials subjected to thermodynamically coupled loading.