The enhanced fracture toughness of advanced functional materials, such as transformation-toughened ceramics and shape memory alloys, stems from the dynamic interplay between solid-state phase transformations and crack propagation. In these systems, stress concentrations near a crack tip trigger a local martensitic transformation, creating a shielding mechanism that dissipates energy and relaxes local stresses. However, computational simulation of this coupled phenomenon remains a significant challenge. The difficulty lies in the system's inherent nonlinearity, driven by the competition among bulk chemical, elastic, and fracture energies across disparate length scales. This study establishes a thermodynamically consistent phase-field framework to resolve these coupled kinetics. We employ a time-dependent Ginzburg-Landau (TDGL) equation to describe the spatiotemporal evolution of non-conserved structural order parameters, utilising a multi-well Landau polynomial to capture the symmetry-breaking nature of specific crystallographic variants. This is coupled with a variational phase-field approach to fracture, in which the sharp crack topology is regularised via a diffuse damage variable. To address the nonlinearity and strong coupling, we solve the governing equations using an implicit numerical scheme. We present numerical results focusing on the stress-induced evolution of martensitic microstructures. The simulations demonstrate the system’s ability to minimise energy by forming self-accommodating variant patterns and twin boundaries under complex loading conditions. Furthermore, we outline the theoretical extension of this framework to include damage evolution. Finally, we discuss the computational challenges associated with this approach, specifically addressing the convergence difficulties arising from the interaction between transformation strains and the crack-driving force.