Inclusions trapped in diamonds provide pivotal information for investigating the Earth's interior conditions at the time of diamond formation. The depth at which these inclusions are encapsulated is crucial, shedding light on the intricate history of diamond growth and subsequent exhumation to the Earth’s surface. The entrapment pressures of these inclusions, obtained through elastic geothermobarometry methods (Kohn et al., 2023; Rustioni et al., 2015) reveal the pressures at which they were enclosed within the diamond. However, applying this method requires certain assumptions, notably assuming all volume changes post-entrapment to be elastic. Viscous deformation or diamond cracking subsequent to growth can release inclusion stress, impacting measured pressures which ultimately can lead to underestimation of growth depths(Angel et al., 2022; Rustioni et al., 2015). Due to the fast exhumation of diamonds to the Earth’ surface, viscous relaxation is often assumed to be negligible. Therefore, we focused on the modeling of brittle failures developed in diamonds around their inclusions, in order to explore the conditions at which fractures may occur in diamonds and to evaluate the associated stress release. We utilized Phase-Field Modeling to analyze fracture propagation and quantify pressure drop resulting from brittle fracture. Our implementation involved an ABAQUS UEL equipped with the BFGS quasi-Newtonian monolithic algorithm, utilizing the history field irreversibility approach and AT2 damage model. Our study culminated in a comparative analysis between Phase-Field Modeling and cohesive zone modeling (CZM)-based discrete models (XFEM). Based on our numerical tests, we found that brittle fracture relaxation accounts for less than 9% of total elastic relaxation, indicating its limited significance. Additionally, we explored how inclusion size, shape, fracture strength, and toughness affect fracture initiation and propagation. This study pioneers phase-field modeling in mineral physics, to predict microscale fractures within a 3D multiaxial loading structure. Furthermore, our comparison with XFEM-CZM adds insights into their applicability in inclusion-matrix problems, highlighting existing discrepancies for further discussion as outlined by Wu et al., 2020 .